|
|
Re: mex function for beginners
Posted:
May 5, 2011 5:11 PM
|
|
"Kaustubha Govind" <kgovind@mathworks.com> wrote in message <ipuua1$gp$1@newscl01ah.mathworks.com>...
> "Ahmed Hakim" <general_ahmed@hotmail.com> wrote in message
> news:iptgvk$g3c$1@fred.mathworks.com...
> > Hello I have problems with using MEX function
> > I purchased some C++ files and I want to use them in MATLAB
> >
> > these files are:
> > SAT_Filter.cpp
> > SAT_DE.cpp SAT_Filter.cpp SAT_Kepler.cpp
> > SAT_RefSys.cpp SAT_Time.cpp SAT_VecMat.cpp RTOD.cpp
> >
> > and their header files:
> > SAT_Force.h SAT_DE.h SAT_Filter.h SAT_Kepler.h
> > SAT_RefSys.h SAT_Time.h SAT_VecMat.h Sat_const.h
> >
> > when using command mex -v SAT_Filter.cpp SAT_DE.cpp SAT_Filter.cpp
> > SAT_Kepler.cpp SAT_RefSys.cpp SAT_Time.cpp SAT_VecMat.cpp RTOD.cpp
> >
> > I have these errors:
> >
> >
> > --> link /out:"SAT_Filter.mexw64" /dll /export:mexFunction
> > /LIBPATH:"C:\PROGRA~1\MATLAB\R2009B\extern\lib\win64\microsoft" libmx.lib
> > libmex.lib libmat.lib /MACHINE:X64 kernel32.lib user32.lib gdi32.lib
> > winspool.lib comdlg32.lib advapi32.lib shell32.lib ole32.lib oleaut32.lib
> > uuid.lib odbc32.lib odbccp32.lib /nologo /incremental:NO
> > /implib:"C:\USERS\AHMEDH~1\APPDATA\LOCAL\TEMP\MEX_PY~1\templib.x"
> > /MAP:"SAT_Filter.mexw64.map"
> > @C:\USERS\AHMEDH~1\APPDATA\LOCAL\TEMP\MEX_PY~1\MEX_TMP.RSP
> > C:\USERS\AHMEDH~1\APPDATA\LOCAL\TEMP\MEX_PY~1\SAT_Filter.obj : warning
> > LNK4042: object specified more than once; extras ignored Creating library
> > C:\USERS\AHMEDH~1\APPDATA\LOCAL\TEMP\MEX_PY~1\templib.x and object
> > C:\USERS\AHMEDH~1\APPDATA\LOCAL\TEMP\MEX_PY~1\templib.exp RTOD.obj : error
> > LNK2019: unresolved external symbol "class Vector __cdecl
> > AccelHarmonic(class Vector const &,class Matrix const
> > &,double,double,class Matrix const &,int,int)"
> > (?AccelHarmonic@@YA?AVVector@@AEBV1@AEBVMatrix@@NN1HH@Z) referenced in
> > function "class Vector __cdecl Accel(double,class Vector const &,int)"
> > (?Accel@@YA?AVVector@@NAEBV1@H@Z) SAT_Filter.mexw64 : fatal error LNK1120:
> > 1 unresolved externals
> > C:\PROGRA~1\MATLAB\R2009B\BIN\MEX.PL: Error: Link of 'SAT_Filter.mexw64'
> > failed.
> > ??? Error using ==> mex at 221
> > Unable to complete successfully.
> >
> > So Can anybody help me?
> >
>
> Hi Ahmed,
>
> Are there any library files that you are failing to link against? Where is
> AccelHarmonic defined?
>
> Kaustubha Govind
Hello Kaustubha
AccelHarmonic is defined in SAT_Force.cpp
here is the C++ code for this file
// SAT_Force.cpp
#include <cmath>
#include <math.h>
#include "SAT_Const.h"
#include "SAT_Force.h"
#include "SAT_RefSys.h"
#include "SAT_VecMat.h"
using namespace std; // Workaround for MS C++ V5.0;
// <cmath> includes <math.h> directly without
// embedding it into the "std" namespace
//Local funtions
namespace ll0
{
// Fractional part of a number (y=x-[x])
double Frac (double x) { return x-floor(x); };
}
//------------------------------------------------------------------------------
//
// Sun
//
// Purpose:
//
// Computes the Sun's geocentric position using a low precision
// analytical series
//
// Input/Output:
//
// Mjd_TT Terrestrial Time (Modified Julian Date)
// <return> Solar position vector [m] with respect to the
// mean equator and equinox of J2000 (EME2000, ICRF)
//
//------------------------------------------------------------------------------
Vector Sun (double Mjd_TT)
{
// Constants
const double eps = 23.43929111*Rad; // Obliquity of J2000 ecliptic
const double T = (Mjd_TT-MJD_J2000)/36525.0; // Julian cent. since J2000
// Variables
double L, M, r;
Vector r_Sun(3);
// Mean anomaly, ecliptic longitude and radius
M = pi2 * ll0:: Frac ( 0.9931267 + 99.9973583*T); // [rad]
L = pi2 * ll0::Frac ( 0.7859444 + M/pi2 +
(6892.0*sin(M)+72.0*sin(2.0*M)) / 1296.0e3); // [rad]
r = 149.619e9 - 2.499e9*cos(M) - 0.021e9*cos(2*M); // [m]
// Equatorial position vector
r_Sun = R_x(-eps) * Vector(r*cos(L),r*sin(L),0.0);
return r_Sun;
}
//------------------------------------------------------------------------------
//
// Moon
//
// Purpose:
//
// Computes the Moon's geocentric position using a low precision
// analytical series
//
// Input/Output:
//
// Mjd_TT Terrestrial Time (Modified Julian Date)
// <return> Lunar position vector [m] with respect to the
// mean equator and equinox of J2000 (EME2000, ICRF)
//
//------------------------------------------------------------------------------
Vector Moon (double Mjd_TT)
{
// Constants
const double eps = 23.43929111*Rad; // Obliquity of J2000 ecliptic
const double T = (Mjd_TT-MJD_J2000)/36525.0; // Julian cent. since J2000
// Variables
double L_0, l,lp, F, D, dL, S, h, N;
double L, B, R, cosB;
Vector r_Moon(3);
// Mean elements of lunar orbit
L_0 = ll0::Frac ( 0.606433 + 1336.851344*T ); // Mean longitude [rev]
// w.r.t. J2000 equinox
l = pi2*ll0::Frac ( 0.374897 + 1325.552410*T ); // Moon's mean anomaly [rad]
lp = pi2*ll0::Frac ( 0.993133 + 99.997361*T ); // Sun's mean anomaly [rad]
D = pi2*ll0::Frac ( 0.827361 + 1236.853086*T ); // Diff. long. Moon-Sun [rad]
F = pi2*ll0::Frac ( 0.259086 + 1342.227825*T ); // Argument of latitude
// Ecliptic longitude (w.r.t. equinox of J2000)
dL = +22640*sin(l) - 4586*sin(l-2*D) + 2370*sin(2*D) + 769*sin(2*l)
-668*sin(lp) - 412*sin(2*F) - 212*sin(2*l-2*D) - 206*sin(l+lp-2*D)
+192*sin(l+2*D) - 165*sin(lp-2*D) - 125*sin(D) - 110*sin(l+lp)
+148*sin(l-lp) - 55*sin(2*F-2*D);
L = pi2 * ll0::Frac( L_0 + dL/1296.0e3 ); // [rad]
// Ecliptic latitude
S = F + (dL+412*sin(2*F)+541*sin(lp)) / Arcs;
h = F-2*D;
N = -526*sin(h) + 44*sin(l+h) - 31*sin(-l+h) - 23*sin(lp+h)
+11*sin(-lp+h) - 25*sin(-2*l+F) + 21*sin(-l+F);
B = ( 18520.0*sin(S) + N ) / Arcs; // [rad]
cosB = cos(B);
// Distance [m]
R = 385000e3 - 20905e3*cos(l) - 3699e3*cos(2*D-l) - 2956e3*cos(2*D)
-570e3*cos(2*l) + 246e3*cos(2*l-2*D) - 205e3*cos(lp-2*D)
-171e3*cos(l+2*D) - 152e3*cos(l+lp-2*D);
// Equatorial coordinates
r_Moon = R_x(-eps) * Vector ( R*cos(L)*cosB, R*sin(L)*cosB, R*sin(B) );
return r_Moon;
}
//------------------------------------------------------------------------------
//
// Illumination
//
// Purpose:
//
// Computes the fractional illumination of a spacecraft in the
// vicinity of the Earth assuming a cylindrical shadow model
//
// Input/output:
//
// r Spacecraft position vector [m]
// r_Sun Sun position vector [m]
// <return> Illumination factor:
// nu=0 Spacecraft in Earth shadow
// nu=1 Spacecraft fully illuminated by the Sun
//
//------------------------------------------------------------------------------
double Illumination ( const Vector& r, const Vector& r_Sun )
{
Vector e_Sun = r_Sun / Norm(r_Sun); // Sun direction unit vector
double s = Dot ( r, e_Sun ); // Projection of s/c position
return ( ( s>0 || Norm(r-s*e_Sun)>R_Earth ) ? 1.0 : 0.0 );
}
//------------------------------------------------------------------------------
//
// AccelHarmonic
//
// Purpose:
//
// Computes the acceleration due to the harmonic gravity field of the
// central body
//
// Input/Output:
//
// r Satellite position vector in the inertial system
// E Transformation matrix to body-fixed system
// GM Gravitational coefficient
// R_ref Reference radius
// CS Spherical harmonics coefficients (un-normalized)
// n_max Maximum degree
// m_max Maximum order (m_max<=n_max; m_max=0 for zonals, only)
// <return> Acceleration (a=d^2r/dt^2)
//
//------------------------------------------------------------------------------
Vector AccelHarmonic (const Vector& r, const Matrix& E,
double GM, double R_ref, const Matrix& CS,
int n_max, int m_max )
{
// Local variables
int n,m; // Loop counters
double r_sqr, rho, Fac; // Auxiliary quantities
double x0,y0,z0; // Normalized coordinates
double ax,ay,az; // Acceleration vector
double C,S; // Gravitational coefficients
Vector r_bf(3); // Body-fixed position
Vector a_bf(3); // Body-fixed acceleration
Matrix V(n_max+2,n_max+2); // Harmonic functions
Matrix W(n_max+2,n_max+2); // work array (0..n_max+1,0..n_max+1)
// Body-fixed position
r_bf = E * r;
// Auxiliary quantities
r_sqr = Dot(r_bf,r_bf); // Square of distance
rho = R_ref*R_ref / r_sqr; /// R^2 / r^2
/// These contains R/r^2
x0 = R_ref * r_bf(0) / r_sqr; // Normalized /// where r_bf(0) is x
y0 = R_ref * r_bf(1) / r_sqr; // coordinates /// where r_bf(1) is y
z0 = R_ref * r_bf(2) / r_sqr; /// where r_bf(2) is z
//
// Evaluate harmonic functions
// V_nm = (R_ref/r)^(n+1) * P_nm(sin(phi)) * cos(m*lambda)
// and
// W_nm = (R_ref/r)^(n+1) * P_nm(sin(phi)) * sin(m*lambda)
// up to degree and order n_max+1
//
// Calculate zonal terms V(n,0); set W(n,0)=0.0
V(0,0) = R_ref / sqrt(r_sqr); /// Known equation (3.31)
W(0,0) = 0.0; /// Known equation (3.31)
V(1,0) = z0 * V(0,0); /// For getting V (1,0) when n = 1 and m= 0 using equation (3.30)
W(1,0) = 0.0; /// All values foe Wn0 are Zeros
/// ---------------> Equation (3.30) with m = 0
for (n=2; n<=n_max+1; n++) {
V(n,0) = ( (2*n-1) * z0 * V(n-1,0) - (n-1) * rho * V(n-2,0) ) / n;
W(n,0) = 0.0; /// All values foe Wn0 are Zeros
};
/// Calculate tesseral and sectorial terms
for (m=1; m<=m_max+1; m++) {
// Calculate V(m,m) .. V(n_max+1,m) ---------------> Equation (3.29)
V(m,m) = (2*m-1) * ( x0*V(m-1,m-1) - y0*W(m-1,m-1) );
W(m,m) = (2*m-1) * ( x0*W(m-1,m-1) + y0*V(m-1,m-1) );
if (m<=n_max) {
V(m+1,m) = (2*m+1) * z0 * V(m,m);
W(m+1,m) = (2*m+1) * z0 * W(m,m);
};
/// ---------------> Equation (3.30)
for (n=m+2; n<=n_max+1; n++) {
V(n,m) = ( (2*n-1)*z0*V(n-1,m) - (n+m-1)*rho*V(n-2,m) ) / (n-m);
W(n,m) = ( (2*n-1)*z0*W(n-1,m) - (n+m-1)*rho*W(n-2,m) ) / (n-m);
};
};
//
// Calculate accelerations ax,ay,az
// ---------------> Equation(3.33)
ax = ay = az = 0.0;
for (m=0; m<=m_max; m++)
for (n=m; n<=n_max ; n++)
if (m==0) {
C = CS(n,0); // = C_n,0
ax -= C * V(n+1,1);
ay -= C * W(n+1,1);
az -= (n+1)*C * V(n+1,0);
}
else {
C = CS(n,m); // = C_n,m
S = CS(m-1,n); // = S_n,m
Fac = 0.5 * (n-m+1) * (n-m+2);
ax += + 0.5 * ( - C * V(n+1,m+1) - S * W(n+1,m+1) )
+ Fac * ( + C * V(n+1,m-1) + S * W(n+1,m-1) );
ay += + 0.5 * ( - C * W(n+1,m+1) + S * V(n+1,m+1) )
+ Fac * ( - C * W(n+1,m-1) + S * V(n+1,m-1) );
az += (n-m+1) * ( - C * V(n+1,m) - S * W(n+1,m) );
};
// Body-fixed acceleration
a_bf = (GM/(R_ref*R_ref)) * Vector(ax,ay,az);
// Inertial acceleration
return Transp(E)*a_bf;
}
//------------------------------------------------------------------------------
//
// AccelPointMass
//
// Purpose:
//
// Computes the perturbational acceleration due to a point mass
//
// Input/Output:
//
// r Satellite position vector
// s Point mass position vector
// GM Gravitational coefficient of point mass
// <return> Acceleration (a=d^2r/dt^2)
//
//------------------------------------------------------------------------------
Vector AccelPointMass (const Vector& r, const Vector& s, double GM)
{
Vector d(3);
// Relative position vector of satellite w.r.t. point mass
d = r - s;
// Acceleration
return (-GM) * ( d/pow(Norm(d),3) + s/pow(Norm(s),3) );
}
//------------------------------------------------------------------------------
//
// AccelSolrad
//
// Purpose:
//
// Computes the acceleration due to solar radiation pressure assuming
// the spacecraft surface normal to the Sun direction
//
// Input/Output:
//
// r Spacecraft position vector
// r_Sun Sun position vector
// Area Cross-section
// mass Spacecraft mass
// CR Solar radiation pressure coefficient
// P0 Solar radiation pressure at 1 AU
// AU Length of one Astronomical Unit
// <return> Acceleration (a=d^2r/dt^2)
//
// Notes:
//
// r, r_sun, Area, mass, P0 and AU must be given in consistent units,
// e.g. m, m^2, kg and N/m^2.
//
//------------------------------------------------------------------------------
Vector AccelSolrad (const Vector& r, const Vector& r_Sun,
double Area, double mass, double CR,
double P0, double AU )
{
Vector d(3);
// Relative position vector of spacecraft w.r.t. Sun
d = r - r_Sun;
// Acceleration
return CR*(Area/mass)*P0*(AU*AU) * d / pow(Norm(d),3);
}
//------------------------------------------------------------------------------
//
// AccelDrag
//
// Purpose:
//
// Computes the acceleration due to the atmospheric drag.
//
// Input/Output:
//
// Mjd_TT Terrestrial Time (Modified Julian Date)
// r Satellite position vector in the inertial system [m]
// v Satellite velocity vector in the inertial system [m/s]
// T Transformation matrix to true-of-date inertial system
// Area Cross-section [m^2]
// mass Spacecraft mass [kg]
// CD Drag coefficient
// <return> Acceleration (a=d^2r/dt^2) [m/s^2]
//
//------------------------------------------------------------------------------
Vector AccelDrag ( double Mjd_TT, const Vector& r, const Vector& v,
const Matrix& T,
double Area, double mass, double CD )
{
// Constants
// Earth angular velocity vector [rad/s]
const double Data_omega[3]= { 0.0, 0.0, 7.29212e-5 };
const Vector omega ( &Data_omega[0], 3);
// Variables
double v_abs, dens;
Vector r_tod(3), v_tod(3);
Vector v_rel(3), a_tod(3);
Matrix T_trp(3,3);
// Transformation matrix to ICRF/EME2000 system
T_trp = Transp(T);
// Position and velocity in true-of-date system
r_tod = T * r;
v_tod = T * v;
// Velocity relative to the Earth's atmosphere
v_rel = v_tod - Cross(omega,r_tod);
v_abs = Norm(v_rel);
// Atmospheric density due to modified Harris-Priester model
dens = Density_HP (Mjd_TT,r_tod);
// Acceleration
a_tod = -0.5*CD*(Area/mass)*dens*v_abs*v_rel;
return T_trp * a_tod;
}
//------------------------------------------------------------------------------
//
// Density_HP
//
// Purpose:
//
// Computes the atmospheric density for the modified Harris-Priester model.
//
// Input/Output:
//
// Mjd_TT Terrestrial Time (Modified Julian Date)
// r_tod Satellite position vector in the inertial system [m]
// <return> Density [kg/m^3]
//
//------------------------------------------------------------------------------
double Density_HP ( double Mjd_TT, const Vector& r_tod )
{
// Constants
const double upper_limit = 1000.0; // Upper height limit [km]
const double lower_limit = 100.0; // Lower height limit [km]
const double ra_lag = 0.523599; // Right ascension lag [rad]
const int n_prm = 3; // Harris-Priester parameter
// 2(6) low(high) inclination
// Harris-Priester atmospheric density model parameters
// Height [km], minimum density, maximum density [gm/km^3]
const int N_Coef = 50;
const double Data_h[N_Coef]= {
100.0, 120.0, 130.0, 140.0, 150.0, 160.0, 170.0, 180.0, 190.0, 200.0,
210.0, 220.0, 230.0, 240.0, 250.0, 260.0, 270.0, 280.0, 290.0, 300.0,
320.0, 340.0, 360.0, 380.0, 400.0, 420.0, 440.0, 460.0, 480.0, 500.0,
520.0, 540.0, 560.0, 580.0, 600.0, 620.0, 640.0, 660.0, 680.0, 700.0,
720.0, 740.0, 760.0, 780.0, 800.0, 840.0, 880.0, 920.0, 960.0,1000.0};
const double Data_c_min[N_Coef] = {
4.974e+05, 2.490e+04, 8.377e+03, 3.899e+03, 2.122e+03, 1.263e+03,
8.008e+02, 5.283e+02, 3.617e+02, 2.557e+02, 1.839e+02, 1.341e+02,
9.949e+01, 7.488e+01, 5.709e+01, 4.403e+01, 3.430e+01, 2.697e+01,
2.139e+01, 1.708e+01, 1.099e+01, 7.214e+00, 4.824e+00, 3.274e+00,
2.249e+00, 1.558e+00, 1.091e+00, 7.701e-01, 5.474e-01, 3.916e-01,
2.819e-01, 2.042e-01, 1.488e-01, 1.092e-01, 8.070e-02, 6.012e-02,
4.519e-02, 3.430e-02, 2.632e-02, 2.043e-02, 1.607e-02, 1.281e-02,
1.036e-02, 8.496e-03, 7.069e-03, 4.680e-03, 3.200e-03, 2.210e-03,
1.560e-03, 1.150e-03 };
const double Data_c_max[N_Coef] = {
4.974e+05, 2.490e+04, 8.710e+03, 4.059e+03, 2.215e+03, 1.344e+03,
8.758e+02, 6.010e+02, 4.297e+02, 3.162e+02, 2.396e+02, 1.853e+02,
1.455e+02, 1.157e+02, 9.308e+01, 7.555e+01, 6.182e+01, 5.095e+01,
4.226e+01, 3.526e+01, 2.511e+01, 1.819e+01, 1.337e+01, 9.955e+00,
7.492e+00, 5.684e+00, 4.355e+00, 3.362e+00, 2.612e+00, 2.042e+00,
1.605e+00, 1.267e+00, 1.005e+00, 7.997e-01, 6.390e-01, 5.123e-01,
4.121e-01, 3.325e-01, 2.691e-01, 2.185e-01, 1.779e-01, 1.452e-01,
1.190e-01, 9.776e-02, 8.059e-02, 5.741e-02, 4.210e-02, 3.130e-02,
2.360e-02, 1.810e-02 };
const Vector h ( &Data_h[0], N_Coef);
const Vector c_min ( &Data_c_min[0], N_Coef);
const Vector c_max ( &Data_c_max[0], N_Coef);
// Variables
int i, ih; // Height section variables
double height; // Earth flattening
double dec_Sun, ra_Sun, c_dec; // Sun declination, right asc.
double c_psi2; // Harris-Priester modification
double density, h_min, h_max, d_min, d_max;// Height, density parameters
Vector r_Sun(3); // Sun position
Vector u(3); // Apex of diurnal bulge
// Satellite height
height = Geodetic(r_tod).h/1000.0; // [km]
// Exit with zero density outside height model limits
if ( height >= upper_limit || height <= lower_limit )
{ return 0.0;
}
// Sun right ascension, declination
r_Sun = Sun ( Mjd_TT );
ra_Sun = atan2( r_Sun(1), r_Sun(0) );
dec_Sun = atan2( r_Sun(2), sqrt( pow(r_Sun(0),2)+pow(r_Sun(1),2) ) );
// Unit vector u towards the apex of the diurnal bulge
// in inertial geocentric coordinates
c_dec = cos(dec_Sun);
u(0) = c_dec * cos(ra_Sun + ra_lag);
u(1) = c_dec * sin(ra_Sun + ra_lag);
u(2) = sin(dec_Sun);
// Cosine of half angle between satellite position vector and
// apex of diurnal bulge
c_psi2 = 0.5 + 0.5 * Dot(r_tod,u)/Norm(r_tod);
// Height index search and exponential density interpolation
ih = 0; // section index reset
for ( i=0; i<=N_Coef-1; i++) // loop over N_Coef height regimes
{
if ( height >= h(i) && height < h(i+1) )
{
ih = i; // ih identifies height section
break;
}
}
h_min = ( h(ih) - h(ih+1) )/log( c_min(ih+1)/c_min(ih) );
h_max = ( h(ih) - h(ih+1) )/log( c_max(ih+1)/c_max(ih) );
d_min = c_min(ih) * exp( (h(ih)-height)/h_min );
d_max = c_max(ih) * exp( (h(ih)-height)/h_max );
// Density computation
density = d_min + (d_max-d_min)*pow(c_psi2,n_prm);
return density * 1.0e-12; // [kg/m^3]
}
//------------------------------------------------------------------------------
//
// AccelMain
//
// Purpose:
//
// Computes the acceleration of an Earth orbiting satellite due to
// - the Earth's harmonic gravity field,
// - the gravitational perturbations of the Sun and Moon
// - the solar radiation pressure and
// - the atmospheric drag
//
// Input/Output:
//
// Mjd_TT Terrestrial Time (Modified Julian Date)
// r Satellite position vector in the ICRF/EME2000 system
// v Satellite velocity vector in the ICRF/EME2000 system
// Area Cross-section
// mass Spacecraft mass
// CR Radiation pressure coefficient
// CD Drag coefficient
// <return> Acceleration (a=d^2r/dt^2) in the ICRF/EME2000 system
//
//------------------------------------------------------------------------------
Vector AccelMain ( double Mjd_TT, const Vector& r, const Vector& v,
double Area, double mass, double CR, double CD )
{
double Mjd_UT1;
Vector a(3), r_Sun(3), r_Moon(3);
Matrix T(3,3), E(3,3);
// Acceleration due to harmonic gravity field
Mjd_UT1 = Mjd_TT;
T = NutMatrix(Mjd_TT) * PrecMatrix(MJD_J2000,Mjd_TT);
E = GHAMatrix(Mjd_UT1) * T;
a = AccelHarmonic ( r,E, Grav.GM,Grav.R_ref,Grav.CS, Grav.n_max,Grav.m_max );
// Luni-solar perturbations
r_Sun = Sun(Mjd_TT);
r_Moon = Moon(Mjd_TT);
a += AccelPointMass ( r, r_Sun, GM_Sun );
a += AccelPointMass ( r, r_Moon, GM_Moon );
// Solar radiation pressure
a += Illumination ( r, r_Sun )
* AccelSolrad ( r, r_Sun, Area, mass, CR, P_Sol, AU );
// Atmospheric drag
a += AccelDrag ( Mjd_TT, r, v, T, Area, mass, CD );
// Acceleration
return a;
}
Thanks
Hope we can find solution
|
|
