gcc Versionsunterschied -> Problem

Dieses Thema im Forum "C/C++" wurde erstellt von upperlimit, 28.06.2006.

  1. #1 upperlimit, 28.06.2006
    upperlimit

    upperlimit Jungspund

    Dabei seit:
    11.04.2005
    Beiträge:
    16
    Zustimmungen:
    0
    Hallo,

    ich habe folgendes Problem:

    Unter SuSE 10.1 (gcc 4.1.0, 64Bit, Kernel 2.6.16.13-4-default) laesst sich die Datei

    Code:
    //------------------------------------------------------------------------------
    //
    // Exercise_4_2.cpp
    // 
    // Purpose: 
    //
    //   Satellite Orbits - Models, Methods, and Applications
    //   Exercise 4-2: Gauss-Jackson 4th-order predictor 
    //
    // Notes:
    //
    //   This software is protected by national and international copyright. 
    //   Any unauthorized use, reproduction or modificaton is unlawful and 
    //   will be prosecuted. Commercial and non-private application of the 
    //   software in any form is strictly prohibited unless otherwise granted
    //   by the authors.
    //
    //   The code is provided without any warranty; without even the implied 
    //   warranty of merchantibility or fitness for a particular purpose.
    //
    // Last modified:
    //
    //   2000/03/04  OMO  Final version (1st edition)
    //   2001/08/17  OMO  Minor upgrades
    //
    // (c) 1999-2001  O. Montenbruck, E. Gill
    //
    //------------------------------------------------------------------------------
    
    #include <cmath>
    #include <iostream>
    #include <iomanip>
    
    #include "SAT_Const.h"
    #include "SAT_Kepler.h"
    #include "SAT_VecMat.h"
    
    using namespace std;
    
    
    //------------------------------------------------------------------------------
    //
    // GJ4P class 
    //
    //------------------------------------------------------------------------------
    
    // Function prototype for second order differential equations
    // void f (double t, const Vector& r, const Vector& v, Vector& a)
    
    typedef void (*GJ4Pfunct)(
      double        t,     // Independent variable
      const Vector& r,     // Position vector 
      const Vector& v,     // Velocity vector  r'=v
      Vector&       a,     // Acceleration     r''=a=f(t,r,v)
      void*         pAux   // Pointer to auxiliary data used within f
    );
    
    
    // Specification
    
    class GJ4P
    {
      public:
    
        // Constructor
        GJ4P (
          GJ4Pfunct   f_,        // Differential equation
          int         n_eqn_,    // Dimension
          void*       pAux_      // Pointer to auxiliary data
          )
        : n_eqn(n_eqn_), f(f_), pAux(pAux_)
        {};
        
        // Initialization step
        void Init (
          double        t_0,     // Initial value of independent variable
          const Vector& r_0,     // Initial value r_0=r(t_0)
          const Vector& v_0,     // Initial value v_0=dr/dt(t_0)
          double        h_       // Step size
        );
    
        // Integration step
        void Step (         
          double&  t,            // Independent variable; updated by t+h
          Vector&  r,            // Value of r(t); updated by r(t+h)
          Vector&  v             // Value of v(t)=dr/dt(t); updated by v(t+h)
        );
    
      private:
    
        // 4th order Runge-Kutta step
        void RK4 (         
          double&  t,            // Independent variable; updated by t+h
          Vector&  r,            // Value of r(t); updated by r(t+h)
          Vector&  v,            // Value of v(t)=dr/dt(t); updated by v(t+h)
          double   h             // Step size
        );
    
        // Elements
        int         n_eqn;       // Dimension
        GJ4Pfunct   f;           // Differential equation
        double      h;           // Step size
        void*       pAux;        // Pointer to auxiliary data requird by f
        Vector      S2,S1;       // First and second sum of acceleration
        Vector      D[4];        // Backward differences of acceleration at t
        Vector      d[4];        // Backward differences of acceleration at t+h
        Vector      r_p,v_p;     // Predictor
    
    };
    
    
    //
    // 4th order Runge-Kutta step for 2nd order differential equation
    //
    
    void GJ4P::RK4 (double&  t, Vector&  r, Vector&  v, double h )
    {
      Vector v_1, v_2, v_3, v_4;
      Vector a_1, a_2, a_3, a_4;
    
      v_1 = v;              f( t      , r            , v_1, a_1, pAux ); 
      v_2 = v+(h/2.0)*a_1;  f( t+h/2.0, r+(h/2.0)*v_1, v_2, a_2, pAux );
      v_3 = v+(h/2.0)*a_2;  f( t+h/2.0, r+(h/2.0)*v_2, v_3, a_3, pAux ); 
      v_4 = v+h*a_3;        f( t+h    , r+h*v_3      , v_4, a_4, pAux );
      
      t = t + h;
      r = r + (h/6.0)*( v_1 + 2.0*v_2 + 2.0*v_3 + v_4 );
      v = v + (h/6.0)*( a_1 + 2.0*a_2 + 2.0*a_3 + a_4 );
    
    };
    
    
    //
    // Initialization of backwards differences from initial conditions
    //
    
    void GJ4P::Init(double t_0, const Vector& r_0, const Vector& v_0, double h_)
    {
      // Order of method
      
      const int m = 4;
    
      // Coefficients gamma/delta of 1st/2nd order Moulton/Cowell corrector method
    
      const double gc[m+1] = {+1.0, -1/2.0, -1/12.0, -1/24.0, -19/720.0 };
      const double dc[m+2] = {+1.0,   -1.0, +1/12.0,     0.0,  -1/240.0, -1/240.0 };
    
      int       i,j;
      double    t = t_0;
      Vector    r = r_0;
      Vector    v = v_0;
    
      // Save step size  
    
      h = h_;     
    
      // Create table of accelerations at past times t-3h, t-2h, and t-h using
      // RK4 steps
    
      f(t,r,v,D[0],pAux);     // D[i]=a(t-ih)
      for (i=1;i<=m-1;i++) {
        RK4(t,r,v,-h);  f(t,r,v,D[i],pAux);   
      };
    
      // Compute backwards differences
      
      for (i=1;i<=m-1;i++) 
        for (j=m-1;j>=i;j--) D[j] = D[j-1]-D[j];
    
      // Initialize backwards sums using 4th order GJ corrector
    
      S1 = v_0/h;              for (i=1;i<=m  ;i++) S1 -= gc[i]*D[i-1];
      S2 = r_0/(h*h)-dc[1]*S1; for (i=2;i<=m+1;i++) S2 -= dc[i]*D[i-2];
    
    };
    
    
    //
    // Step from t to t+h
    //
    
    void GJ4P::Step (double& t, Vector& r, Vector& v) 
    {
      // Order of method
      
      const int m = 4;  
    
      // Coefficients gamma/delta of 1st/2nd order Bashforth/Stoermr predictor
      
      const double gp[m+1] = {+1.0, +1/2.0, +5/12.0,  +3/8.0, +251/720.0 };
      const double dp[m+2] = {+1.0,    0.0, +1/12.0, +1/12.0,  +19/240.0,  +3/40.0 };
    
      int i;
      
      // 4th order predictor
    
      r_p = dp[0]*S2; for(i=2;i<=m+1;i++) r_p += dp[i]*D[i-2]; r_p = (h*h)*r_p;
      v_p = gp[0]*S1; for(i=1;i<=m  ;i++) v_p += gp[i]*D[i-1]; v_p =     h*v_p;
    
      // Update backwards difference table
    
      f ( t+h, r_p,v_p, d[0], pAux );               // Acceleration at t+h
      for (i=1;i<=m-1;i++) d[i]=d[i-1]-D[i-1];      // New differences at t+h
      for (i=0;i<=m-1;i++) D[i]=d[i];               // Update differences 
      S1 += d[0];  S2 += S1;                        // Update sums
    
      // Update independent variable and solution
    
      t = t + h;
      r = r_p;
      v = v_p;
    
    };
    
    
    //------------------------------------------------------------------------------
    //
    // f_Kep3D
    //
    // Purpose:
    // 
    //   Computes the second time derivative of the position vector for the 
    //   normalized (GM=1) Kepler's problem in three dimensions
    //
    // Note:
    //
    //   pAux is expected to point to an integer variable that will be incremented
    //   by one on each call of f_Kep3D
    //
    //------------------------------------------------------------------------------
    
    void f_Kep3D ( double t, const Vector& r, const Vector& v, 
                   Vector& a, void* pAux )
    {
      
      // Pointer to auxiliary integer variable used as function call counter
      int* pCalls = static_cast<int*>(pAux);
      
      // 2nd order derivative d^2(r)/dt^2 of the position vector
      a = -r/(pow(Norm(r),3));
      
      // Increment function call count
      (*pCalls)++;
    
    };
    
    
    //------------------------------------------------------------------------------
    //
    // Main program
    //
    //------------------------------------------------------------------------------
    
    int main() {
    
      // Constants
      
      const double  GM    = 1.0;                   // Gravitational coefficient
      const double  e     = 0.1;                   // Eccentricity
      const double  t_end = 20.0;                  // End time
      const Vector  Kep(1.0,e,0.0,0.0,0.0,0.0);    // (a,e,i,Omega,omega,M)
      const Vector  y_ref = State(GM,Kep,t_end);   // Reference solution
    
      const int     Steps[] = { 100, 300, 600, 1000, 1500, 2000, 3000, 4000 };
      
      // Variables
      
      int     nCalls;                              // Function call count
      int     iCase;
      double  t,h;                                 // Time and step size
      Vector  y(6);                                // State vector
      Vector  r(3);
      Vector  v(3);
    
      GJ4P    Orbit(f_Kep3D,3,&nCalls);            // Object for integrating the
                                                   // 2nd order diff. equation
                                                   // defined by f_Kep3D using the
                                                   // 4th-order GJ predictor
      
      // Header 
    
      cout << "Exercise 4-2: Gauss-Jackson 4th-order predictor" 
           << endl << endl;
      cout << "  Problem D1 (e=0.1)" << endl << endl;
      cout << "  N_fnc   Accuracy   Digits " << endl;
        
      // Loop over test cases
    
      for (iCase=0; iCase<8; iCase++) {
      
        // Step size
        h = t_end/Steps[iCase];
    
        // Initial values
        t = 0.0;
        r = Vector( 1.0-e, 0.0, 0.0 );
        v = Vector( 0.0, sqrt((1+e)/(1-e)), 0.0 );
        nCalls = 0;
    
        // Integration from t=t to t=t_end
        Orbit.Init ( t, r, v, h );
        for (int i=1; i<=Steps[iCase]; i++)  
          Orbit.Step( t, r, v );
        y = Stack(r,v);
    
        // Output
        cout << fixed  << setw(6) << nCalls
             << scientific << setprecision(3) << setw(13)
             << Norm(y-y_ref) 
             << fixed << setprecision(2) << setw(7)
             << -log10(Norm(y-y_ref)) << endl;
      
      };
    
      return 0;
    
    }
    
    ( http://www.fjtakalab.ynu.ac.jp/fjtaka_lec/VLBI/Source_Astro/Source/Exercise_4_2.cpp )

    mit den Header-Files

    http://www.fjtakalab.ynu.ac.jp/fjtaka_lec/VLBI/Source_Astro/Source/SAT_Const.h
    http://www.fjtakalab.ynu.ac.jp/fjtaka_lec/VLBI/Source_Astro/Source/SAT_Kepler.h
    und
    http://www.fjtakalab.ynu.ac.jp/fjtaka_lec/VLBI/Source_Astro/Source/SAT_VecMat.h

    nicht kompilieren.


    Der Befehl

    Code:
    $gcc -I. -c Exercise_4_2.cpp
    
    bzw.
    Code:
    $g++ -I. -c Exercise_4_2.cpp
    
    bringt die Fehlermeldung:

    Code:
    SAT_VecMat.h:106: error: ‘Matrix’ does not name a type
    SAT_VecMat.h:109: error: expected ‘,’ or ‘...’ before ‘&’ token
    SAT_VecMat.h:109: error: ISO C++ forbids declaration of ‘Matrix’ with no type
    SAT_VecMat.h:109: error: ‘Vector operator*(int)’ must have an argument of class or enumerated type
    SAT_VecMat.h:110: error: expected ‘,’ or ‘...’ before ‘&’ token
    SAT_VecMat.h:110: error: ISO C++ forbids declaration of ‘Matrix’ with no type
    SAT_VecMat.h:113: error: ‘Matrix’ does not name a type
    
    Interessanterweise wird dieselbe Aufgabe mit dem gcc 3.3.3 (SuSE 9.0, 32Bit, Kernel 2.4.x-default) anstandslos erledigt.

    Meine Frage lautet also:
    Wie kann ich es erreichen, dass die C++-Quelle auch unter gcc 4.1.0 compiliert wird ? (Um einen Unterschied zw. 32- und 64Bit-Architektur scheint es nicht zu gehen, da dasselbe Problem auch auf einer 32Bit-Maschine auftritt, allerdings auch mit gcc 4.1.0.)

    Jeder Vorschlag ist willkommen.

    Gruesse,
    upperlimit
     
  2. Anzeige

    Schau dir mal diese Kategorie an. Dort findest du bestimmt etwas.
    Registrieren bzw. einloggen, um diese und auch andere Anzeigen zu deaktivieren
  3. #2 Wolfgang, 28.06.2006
    Wolfgang

    Wolfgang Foren Gott

    Dabei seit:
    24.04.2005
    Beiträge:
    3.978
    Zustimmungen:
    0
    Ort:
    Erfurt
    Hallo
    Wenn der Code nunmal für den gcc 3.3.3 geschrieben wurde, dann musst du ihn entweder selbst bearbeiten - der neueren Version anpassen - , oder eben mit dem erforderlichen Compiler compilieren.
    Dafür gibt es ja in der Regel die Hinweise, welche Version benötigt wird. Steht meist in der README.

    Gruß Wolfgang
     
  4. rikola

    rikola Foren Gott

    Dabei seit:
    23.08.2005
    Beiträge:
    2.133
    Zustimmungen:
    0
    Ich wuerde mal sagen, dass der Code nichts Versionsspezifisches enthaelt. gcc-4.0 ist einfach strenger und enger am Standard. In der Datei SAT_VecMat.h wird die Klasse Matrix innerhalb der Klasse Vector benutzt, ohne sie vorher zu deklarieren. Schreibe einfach ueber Zeile 44, also vor der Deklaration der Klasse Vector, ein 'class Matrix;' hin, dann sollte zumindest diese Fehlermeldung verschwinden. Es wundert mich, dass gcc-3.3 das hat durchgehen lassen.
     
  5. #4 upperlimit, 29.06.2006
    upperlimit

    upperlimit Jungspund

    Dabei seit:
    11.04.2005
    Beiträge:
    16
    Zustimmungen:
    0
    Hallo,

    nach einer Parallelinstallation von gcc 4.1.0 und gcc 3.3.6 (SuSE 10.0, 64Bit, Kernel 2.6.16.13-4-default) habe ich es nun hinbekommen, dass die Datei mit dem gcc 3.3.6 kompiliert werden kann.

    Dein Tipp, rikola, mit dem Einfuegen von 'class Matrix;' in die SAT_VecMat.h fuehrt allerdings zu anderen Fehler in der Datei SAT_RefSys1.cpp

    Code:
    //------------------------------------------------------------------------------
    //
    // SAT_RefSys1.cpp
    // 
    // Purpose: 
    //
    //   Transformations betweeen celestial and terrestrial reference systems
    //
    // Notes:
    //
    //   This software is protected by national and international copyright. 
    //   Any unauthorized use, reproduction or modificaton is unlawful and 
    //   will be prosecuted. Commercial and non-private application of the 
    //   software in any form is strictly prohibited unless otherwise granted
    //   by the authors.
    //
    //   The code is provided without any warranty; without even the implied 
    //   warranty of merchantibility or fitness for a particular purpose.
    //
    // Last modified:
    //
    //   2000/03/04  OMO  Final version (1st edition)
    //   
    // (c) 1999-2000  O. Montenbruck, E. Gill
    //
    //------------------------------------------------------------------------------
    
    #include <cmath>
    
    #ifdef __GNUC__   // GNU C++ adaptation
    #include <float.h>
    #else             // Standard C++ version
    #include <limits>
    #endif
    
    #include "SAT_Const.h"
    #include "SAT_RefSys1.h"
    #include "SAT_VecMat.h"
    
    using namespace std;
    
    
    
    //
    // Local declarations
    //
    
    namespace { 
      // Machine accuracy
      #ifdef __GNUC__   // GNU C++ adaptation
      const double eps_mach = DBL_EPSILON;
      #else             // Standard C++ version
      const double eps_mach = numeric_limits<double>::epsilon();
      #endif
      // Fractional part of a number (y=x-[x])
      double Frac (double x) { return x-floor(x); };
      // x mod y
      double Modulo (double x, double y) { return y*Frac(x/y); }
    }            
    
    
    //------------------------------------------------------------------------------
    //
    // IERS (Class implementation)
    //
    // Purpose:
    //
    //   Management of IERS time and polar motion data
    //  
    //   p.157-169
    //   5.1 Time
    //------------------------------------------------------------------------------
    
    // Class constants
    
    const double IERS::TT_TAI  = +32.184; // TT-TAI time difference [s] p.162, 5.3
    const double IERS::GPS_TAI = -19.0;   // GPS-TAI time difference [s] p.162, 5.4
    
    // Default values of Earth rotation parameters
    
    double IERS::UT1_UTC_ = 0.0;          // UT1-UTC time difference [s]
    double IERS::UTC_TAI_ = 0.0;          // UTC-TAI time difference [s] 
    double IERS::x_pole_  = 0.0;          // Pole coordinate [rad]
    double IERS::y_pole_  = 0.0;          // Pole coordinate [rad]
    
    // Setting of IERS Earth rotation parameters
    // (UT1-UTC [s],UTC-TAI [s], x ["], y ["])
    // cf. Fig.5.4, p.166, Buuletin B of the IERS 
    
    void IERS::Set(double UT1_UTC, double UTC_TAI,
                   double x_pole, double y_pole) 
    {
       UT1_UTC_ = UT1_UTC;
       UTC_TAI_ = UTC_TAI;
       x_pole_  = x_pole/Arcs;
       y_pole_  = y_pole/Arcs;
    };
    
    // UT1-UTC time difference [s]
    
    double IERS::UT1_UTC(double Mjd_UTC) { return UT1_UTC_; };
         
    // UTC-TAI time difference [s]
    
    double IERS::UTC_TAI(double Mjd_UTC) { return UTC_TAI_; };
    
    // TT-UTC time difference [s]
    
    double IERS::TT_UTC(double Mjd_UTC) { return TT_TAI-UTC_TAI_; };
    
    // GPS-UTC time difference [s]
    
    double IERS::GPS_UTC(double Mjd_UTC) { return GPS_TAI-UTC_TAI_; };
    
    // Pole coordinate [rad]
    
    double IERS::x_pole(double Mjd_UTC) { return x_pole_; };
    
    // Pole coordinate [rad]
    
    double IERS::y_pole(double Mjd_UTC) { return y_pole_; };
    
    
    
    //------------------------------------------------------------------------------
    //
    // MeanObliquity
    //
    // Purpose:
    //
    //   Computes the mean obliquity of the ecliptic
    //
    // Input/Output:
    //
    //   Mjd_TT    Modified Julian Date (Terrestrial Time)
    //   <return>  Mean obliquity of the ecliptic
    //
    //  p.176
    //  5.3.2 Coordinate Changes due to Precession
    //  5.42
    //------------------------------------------------------------------------------
    
    double MeanObliquity (double Mjd_TT)
    {
      const double T = (Mjd_TT-MJD_J2000)/36525.0;
    
      return 
        Rad *( 23.43929111-(46.8150+(0.00059-0.001813*T)*T)*T/3600.0 );
    }
    
    
    //------------------------------------------------------------------------------
    //
    // EclMatrix
    //
    // Purpose:
    //
    //   Transformation of equatorial to ecliptical coordinates
    //
    // Input/Output:
    //
    //   Mjd_TT    Modified Julian Date (Terrestrial Time)
    //   <return>  Transformation matrix
    //
    //  p.170
    //  5.2 Celestial and Terrestrial Reference Systems
    //  5.21
    //------------------------------------------------------------------------------
    
    Matrix EclMatrix (double Mjd_TT)
    {
      return R_x ( MeanObliquity(Mjd_TT) );
    }
    
    
    //------------------------------------------------------------------------------
    //
    // PrecMatrix
    //
    // Purpose:
    //
    //   Precession transformation of equatorial coordinates
    //
    // Input/Output:
    //
    //   Mjd_1     Epoch given (Modified Julian Date TT)
    //   MjD_2     Epoch to precess to (Modified Julian Date TT)
    //   <return>  Precession transformation matrix
    //
    //  p.174-177
    //  5.3.2 Coordinate Changes due to Precession
    //  5.44,5.46,5.47
    //------------------------------------------------------------------------------
    
    Matrix PrecMatrix (double Mjd_1, double Mjd_2)
    {
    
      // Constants
    
      const double T  = (Mjd_1-MJD_J2000)/36525.0;
      const double dT = (Mjd_2-Mjd_1)/36525.0;
      
      // Variables
    
      double zeta,z,theta;
    
      // Precession angles
      
      zeta  =  ( (2306.2181+(1.39656-0.000139*T)*T)+
                    ((0.30188-0.000344*T)+0.017998*dT)*dT )*dT/Arcs;
      z     =  zeta + ( (0.79280+0.000411*T)+0.000205*dT)*dT*dT/Arcs;
      theta =  ( (2004.3109-(0.85330+0.000217*T)*T)-
                    ((0.42665+0.000217*T)+0.041833*dT)*dT )*dT/Arcs;
    
      // Precession matrix
      
      return R_z(-z) * R_y(theta) * R_z(-zeta);
    
    }    
    
    
    //------------------------------------------------------------------------------
    //
    // NutAngles 
    //
    // Purpose:
    //
    //   Nutation in longitude and obliquity
    //
    // Input/Output:
    //
    //   Mjd_TT    Modified Julian Date (Terrestrial Time)
    //   <return>  Nutation matrix
    //
    //  p.178-181
    //  5.3.3 Nutation
    //  Table 5.2, p.179, the IAU 1980 nutation theory, 5.58,5.59,5.60
    //------------------------------------------------------------------------------
    
    void NutAngles (double Mjd_TT, double& dpsi, double& deps)
    {
    
      // Constants
    
      const double T  = (Mjd_TT-MJD_J2000)/36525.0;
      const double T2 = T*T;
      const double T3 = T2*T;
      const double rev = 360.0*3600.0;  // arcsec/revolution
    
      const int  N_coeff = 106;
      const long C[N_coeff][9] =
      {
       //
       // l  l' F  D Om    dpsi    *T     deps     *T       #
       //
        {  0, 0, 0, 0, 1,-1719960,-1742,  920250,   89 },   //   1
        {  0, 0, 0, 0, 2,   20620,    2,   -8950,    5 },   //   2
        { -2, 0, 2, 0, 1,     460,    0,    -240,    0 },   //   3
        {  2, 0,-2, 0, 0,     110,    0,       0,    0 },   //   4
        { -2, 0, 2, 0, 2,     -30,    0,      10,    0 },   //   5
        {  1,-1, 0,-1, 0,     -30,    0,       0,    0 },   //   6
        {  0,-2, 2,-2, 1,     -20,    0,      10,    0 },   //   7
        {  2, 0,-2, 0, 1,      10,    0,       0,    0 },   //   8
        {  0, 0, 2,-2, 2, -131870,  -16,   57360,  -31 },   //   9
        {  0, 1, 0, 0, 0,   14260,  -34,     540,   -1 },   //  10
        {  0, 1, 2,-2, 2,   -5170,   12,    2240,   -6 },   //  11
        {  0,-1, 2,-2, 2,    2170,   -5,    -950,    3 },   //  12
        {  0, 0, 2,-2, 1,    1290,    1,    -700,    0 },   //  13
        {  2, 0, 0,-2, 0,     480,    0,      10,    0 },   //  14
        {  0, 0, 2,-2, 0,    -220,    0,       0,    0 },   //  15
        {  0, 2, 0, 0, 0,     170,   -1,       0,    0 },   //  16
        {  0, 1, 0, 0, 1,    -150,    0,      90,    0 },   //  17
        {  0, 2, 2,-2, 2,    -160,    1,      70,    0 },   //  18
        {  0,-1, 0, 0, 1,    -120,    0,      60,    0 },   //  19
        { -2, 0, 0, 2, 1,     -60,    0,      30,    0 },   //  20
        {  0,-1, 2,-2, 1,     -50,    0,      30,    0 },   //  21
        {  2, 0, 0,-2, 1,      40,    0,     -20,    0 },   //  22
        {  0, 1, 2,-2, 1,      40,    0,     -20,    0 },   //  23
        {  1, 0, 0,-1, 0,     -40,    0,       0,    0 },   //  24
        {  2, 1, 0,-2, 0,      10,    0,       0,    0 },   //  25
        {  0, 0,-2, 2, 1,      10,    0,       0,    0 },   //  26
        {  0, 1,-2, 2, 0,     -10,    0,       0,    0 },   //  27
        {  0, 1, 0, 0, 2,      10,    0,       0,    0 },   //  28
        { -1, 0, 0, 1, 1,      10,    0,       0,    0 },   //  29
        {  0, 1, 2,-2, 0,     -10,    0,       0,    0 },   //  30
        {  0, 0, 2, 0, 2,  -22740,   -2,    9770,   -5 },   //  31
        {  1, 0, 0, 0, 0,    7120,    1,     -70,    0 },   //  32
        {  0, 0, 2, 0, 1,   -3860,   -4,    2000,    0 },   //  33
        {  1, 0, 2, 0, 2,   -3010,    0,    1290,   -1 },   //  34
        {  1, 0, 0,-2, 0,   -1580,    0,     -10,    0 },   //  35
        { -1, 0, 2, 0, 2,    1230,    0,    -530,    0 },   //  36
        {  0, 0, 0, 2, 0,     630,    0,     -20,    0 },   //  37
        {  1, 0, 0, 0, 1,     630,    1,    -330,    0 },   //  38
        { -1, 0, 0, 0, 1,    -580,   -1,     320,    0 },   //  39
        { -1, 0, 2, 2, 2,    -590,    0,     260,    0 },   //  40
        {  1, 0, 2, 0, 1,    -510,    0,     270,    0 },   //  41
        {  0, 0, 2, 2, 2,    -380,    0,     160,    0 },   //  42
        {  2, 0, 0, 0, 0,     290,    0,     -10,    0 },   //  43
        {  1, 0, 2,-2, 2,     290,    0,    -120,    0 },   //  44
        {  2, 0, 2, 0, 2,    -310,    0,     130,    0 },   //  45
        {  0, 0, 2, 0, 0,     260,    0,     -10,    0 },   //  46
        { -1, 0, 2, 0, 1,     210,    0,    -100,    0 },   //  47
        { -1, 0, 0, 2, 1,     160,    0,     -80,    0 },   //  48
        {  1, 0, 0,-2, 1,    -130,    0,      70,    0 },   //  49
        { -1, 0, 2, 2, 1,    -100,    0,      50,    0 },   //  50
        {  1, 1, 0,-2, 0,     -70,    0,       0,    0 },   //  51
        {  0, 1, 2, 0, 2,      70,    0,     -30,    0 },   //  52
        {  0,-1, 2, 0, 2,     -70,    0,      30,    0 },   //  53
        {  1, 0, 2, 2, 2,     -80,    0,      30,    0 },   //  54
        {  1, 0, 0, 2, 0,      60,    0,       0,    0 },   //  55
        {  2, 0, 2,-2, 2,      60,    0,     -30,    0 },   //  56
        {  0, 0, 0, 2, 1,     -60,    0,      30,    0 },   //  57
        {  0, 0, 2, 2, 1,     -70,    0,      30,    0 },   //  58
        {  1, 0, 2,-2, 1,      60,    0,     -30,    0 },   //  59
        {  0, 0, 0,-2, 1,     -50,    0,      30,    0 },   //  60
        {  1,-1, 0, 0, 0,      50,    0,       0,    0 },   //  61
        {  2, 0, 2, 0, 1,     -50,    0,      30,    0 },   //  62
        {  0, 1, 0,-2, 0,     -40,    0,       0,    0 },   //  63
        {  1, 0,-2, 0, 0,      40,    0,       0,    0 },   //  64
        {  0, 0, 0, 1, 0,     -40,    0,       0,    0 },   //  65
        {  1, 1, 0, 0, 0,     -30,    0,       0,    0 },   //  66
        {  1, 0, 2, 0, 0,      30,    0,       0,    0 },   //  67
        {  1,-1, 2, 0, 2,     -30,    0,      10,    0 },   //  68
        { -1,-1, 2, 2, 2,     -30,    0,      10,    0 },   //  69
        { -2, 0, 0, 0, 1,     -20,    0,      10,    0 },   //  70
        {  3, 0, 2, 0, 2,     -30,    0,      10,    0 },   //  71
        {  0,-1, 2, 2, 2,     -30,    0,      10,    0 },   //  72
        {  1, 1, 2, 0, 2,      20,    0,     -10,    0 },   //  73
        { -1, 0, 2,-2, 1,     -20,    0,      10,    0 },   //  74
        {  2, 0, 0, 0, 1,      20,    0,     -10,    0 },   //  75
        {  1, 0, 0, 0, 2,     -20,    0,      10,    0 },   //  76
        {  3, 0, 0, 0, 0,      20,    0,       0,    0 },   //  77
        {  0, 0, 2, 1, 2,      20,    0,     -10,    0 },   //  78
        { -1, 0, 0, 0, 2,      10,    0,     -10,    0 },   //  79
        {  1, 0, 0,-4, 0,     -10,    0,       0,    0 },   //  80
        { -2, 0, 2, 2, 2,      10,    0,     -10,    0 },   //  81
        { -1, 0, 2, 4, 2,     -20,    0,      10,    0 },   //  82
        {  2, 0, 0,-4, 0,     -10,    0,       0,    0 },   //  83
        {  1, 1, 2,-2, 2,      10,    0,     -10,    0 },   //  84
        {  1, 0, 2, 2, 1,     -10,    0,      10,    0 },   //  85
        { -2, 0, 2, 4, 2,     -10,    0,      10,    0 },   //  86
        { -1, 0, 4, 0, 2,      10,    0,       0,    0 },   //  87
        {  1,-1, 0,-2, 0,      10,    0,       0,    0 },   //  88
        {  2, 0, 2,-2, 1,      10,    0,     -10,    0 },   //  89
        {  2, 0, 2, 2, 2,     -10,    0,       0,    0 },   //  90
        {  1, 0, 0, 2, 1,     -10,    0,       0,    0 },   //  91
        {  0, 0, 4,-2, 2,      10,    0,       0,    0 },   //  92
        {  3, 0, 2,-2, 2,      10,    0,       0,    0 },   //  93
        {  1, 0, 2,-2, 0,     -10,    0,       0,    0 },   //  94
        {  0, 1, 2, 0, 1,      10,    0,       0,    0 },   //  95
        { -1,-1, 0, 2, 1,      10,    0,       0,    0 },   //  96
        {  0, 0,-2, 0, 1,     -10,    0,       0,    0 },   //  97
        {  0, 0, 2,-1, 2,     -10,    0,       0,    0 },   //  98
        {  0, 1, 0, 2, 0,     -10,    0,       0,    0 },   //  99
        {  1, 0,-2,-2, 0,     -10,    0,       0,    0 },   // 100
        {  0,-1, 2, 0, 1,     -10,    0,       0,    0 },   // 101
        {  1, 1, 0,-2, 1,     -10,    0,       0,    0 },   // 102
        {  1, 0,-2, 2, 0,     -10,    0,       0,    0 },   // 103
        {  2, 0, 0, 2, 0,      10,    0,       0,    0 },   // 104
        {  0, 0, 2, 4, 2,     -10,    0,       0,    0 },   // 105
        {  0, 1, 0, 1, 0,      10,    0,       0,    0 }    // 106
       };
    
      // Variables
    
      double  l, lp, F, D, Om;
      double  arg;
      
    
      // Mean arguments of luni-solar motion
      //
      //   l   mean anomaly of the Moon
      //   l'  mean anomaly of the Sun
      //   F   mean argument of latitude
      //   D   mean longitude elongation of the Moon from the Sun 
      //   Om  mean longitude of the ascending node
      
      l  = Modulo (  485866.733 + (1325.0*rev +  715922.633)*T    
                                     + 31.310*T2 + 0.064*T3, rev );
      lp = Modulo ( 1287099.804 + (  99.0*rev + 1292581.224)*T             
                                     -  0.577*T2 - 0.012*T3, rev );
      F  = Modulo (  335778.877 + (1342.0*rev +  295263.137)*T    
                                     - 13.257*T2 + 0.011*T3, rev );
      D  = Modulo ( 1072261.307 + (1236.0*rev + 1105601.328)*T    
                                     -  6.891*T2 + 0.019*T3, rev );
      Om = Modulo (  450160.280 - (   5.0*rev +  482890.539)*T    
                                     +  7.455*T2 + 0.008*T3, rev );
    
      // Nutation in longitude and obliquity [rad]
    
      deps = dpsi = 0.0;
      for (int i=0; i<N_coeff; i++) {
        arg  =  ( C[i][0]*l+C[i][1]*lp+C[i][2]*F+C[i][3]*D+C[i][4]*Om ) / Arcs;
        dpsi += ( C[i][5]+C[i][6]*T ) * sin(arg);
        deps += ( C[i][7]+C[i][8]*T ) * cos(arg);
      };
          
      dpsi = 1.0E-5 * dpsi/Arcs;
      deps = 1.0E-5 * deps/Arcs;
    
    }
    
    
    //------------------------------------------------------------------------------
    //
    // NutMatrix 
    //
    // Purpose:
    //
    //   Transformation from mean to true equator and equinox
    //
    // Input/Output:
    //
    //   Mjd_TT    Modified Julian Date (Terrestrial Time)
    //   <return>  Nutation matrix
    //
    //  p.178-181
    //  5.3.3 Nutation
    //  5.62
    //------------------------------------------------------------------------------
    
    Matrix NutMatrix (double Mjd_TT)
    {
    
      double dpsi, deps, eps;
      
      // Mean obliquity of the ecliptic
    
      eps = MeanObliquity(Mjd_TT);
    
      // Nutation in longitude and obliquity
    
      NutAngles (Mjd_TT, dpsi,deps);
    
      // Transformation from mean to true equator and equinox
    
      return  R_x(-eps-deps)*R_z(-dpsi)*R_x(+eps);
    
    }
    
    
    //------------------------------------------------------------------------------
    //
    // NutMatrixSimple
    //
    // Purpose:
    //
    //   Transformation from mean to true equator and equinox (low precision)
    //
    // Input/Output:
    //
    //   Mjd_TT    Modified Julian Date (Terrestrial Time)
    //   <return>  Nutation matrix
    //
    //------------------------------------------------------------------------------
    
    Matrix NutMatrixSimple (double Mjd_TT)
    {
    
      // Constants
    
      const double T  = (Mjd_TT-MJD_J2000)/36525.0;
    
      // Variables
    
      double  ls, D, F, N;
      double  eps, dpsi, deps;
    
      // Mean arguments of luni-solar motion
      
      ls = pi2*Frac(0.993133+  99.997306*T);   // mean anomaly Sun          
      D  = pi2*Frac(0.827362+1236.853087*T);   // diff. longitude Moon-Sun  
      F  = pi2*Frac(0.259089+1342.227826*T);   // mean argument of latitude 
      N  = pi2*Frac(0.347346-   5.372447*T);   // longit. ascending node    
    
      // Nutation angles
      
      dpsi = ( -17.200*sin(N)   - 1.319*sin(2*(F-D+N)) - 0.227*sin(2*(F+N))
               + 0.206*sin(2*N) + 0.143*sin(ls) ) / Arcs;
      deps = ( + 9.203*cos(N)   + 0.574*cos(2*(F-D+N)) + 0.098*cos(2*(F+N))
               - 0.090*cos(2*N)                 ) / Arcs;
    
      // Mean obliquity of the ecliptic
    
      eps  = 0.4090928-2.2696E-4*T;    
    
      return  R_x(-eps-deps)*R_z(-dpsi)*R_x(+eps);
    
    }
    
    
    //------------------------------------------------------------------------------
    //
    // EqnEquinox 
    //
    // Purpose:
    //
    //   Computation of the equation of the equinoxes
    //
    // Input/Output:
    //
    //   Mjd_TT    Modified Julian Date (Terrestrial Time)
    //   <return>  Equation of the equinoxes
    //
    // Notes:
    //
    //   The equation of the equinoxes dpsi*cos(eps) is the right ascension of the 
    //   mean equinox referred to the true equator and equinox and is equal to the 
    //   difference between apparent and mean sidereal time.
    //
    //  p.181
    //  5.4.1 Rotation About the Celestial Ephemeris Pole
    //  5.66
    //------------------------------------------------------------------------------
    
    double EqnEquinox (double Mjd_TT)
    {
      double dpsi, deps;              // Nutation angles
    
      // Nutation in longitude and obliquity 
          
      NutAngles (Mjd_TT, dpsi,deps );
    
      // Equation of the equinoxes
    
      return  dpsi * cos ( MeanObliquity(Mjd_TT) );
    
    };
    
    
    //------------------------------------------------------------------------------
    //
    // GMST
    //
    // Purpose:
    //
    //   Greenwich Mean Sidereal Time
    //
    // Input/Output:
    //
    //   Mjd_UT1   Modified Julian Date UT1
    //   <return>  GMST in [rad]
    //
    //  p.165-169
    //  Sidereal Time and Universal Time
    //  5.15 -5.20
    //------------------------------------------------------------------------------
    
    double GMST (double Mjd_UT1)
    {
    
      // Constants
    
      const double Secs = 86400.0;        // Seconds per day
    
      // Variables
    
      double Mjd_0,UT1,T_0,T,gmst;
    
      // Mean Sidereal Time
      
      Mjd_0 = floor(Mjd_UT1);
      UT1   = Secs*(Mjd_UT1-Mjd_0);          // [s]
      T_0   = (Mjd_0  -MJD_J2000)/36525.0; 
      T     = (Mjd_UT1-MJD_J2000)/36525.0; 
    
      gmst  = 24110.54841 + 8640184.812866*T_0 + 1.002737909350795*UT1
              + (0.093104-6.2e-6*T)*T*T; // [s]
    
      return  pi2*Frac(gmst/Secs);       // [rad], 0..2pi
    
    }
    
    
    //------------------------------------------------------------------------------
    //
    // GAST
    //
    // Purpose:
    //
    //   Greenwich Apparent Sidereal Time
    //
    // Input/Output:
    //
    //   Mjd_UT1   Modified Julian Date UT1
    //   <return>  GMST in [rad]
    //
    //  p.181
    //  5.4.1 Rotation About the Celestial Ephemeris Pole
    //  5.66
    //------------------------------------------------------------------------------
    
    double GAST (double Mjd_UT1)
    {
      return Modulo ( GMST(Mjd_UT1) + EqnEquinox(Mjd_UT1), pi2 );
    }
    
    
    //------------------------------------------------------------------------------
    //
    // GHAMatrix
    //
    // Purpose:
    //
    //   Transformation from true equator and equinox to Earth equator and 
    //   Greenwich meridian system 
    //
    // Input/Output:
    //
    //   Mjd_UT1   Modified Julian Date UT1
    //   <return>  Greenwich Hour Angle matrix
    //
    //  p.181
    //  5.4.1 Rotation About the Celestial Ephemeris Pole
    //  5.67
    //------------------------------------------------------------------------------
    
    Matrix GHAMatrix (double Mjd_UT1)
    {
      return  R_z ( GAST(Mjd_UT1) );
    }
    
    
    //------------------------------------------------------------------------------
    //
    // PoleMatrix
    //
    // Purpose:
    //
    //   Transformation from pseudo Earth-fixed to Earth-fixed coordinates
    //   for a given date
    //
    // Input/Output:
    //
    //   Mjd_UTC   Modified Julian Date UTC
    //   <return>  Pole matrix
    //
    //  p.185
    //  5.4.4 Transformation to the International Reference Pole
    //  5.75
    //------------------------------------------------------------------------------
    
    Matrix PoleMatrix (double Mjd_UTC)
    {
       return  R_y(-IERS::x_pole(Mjd_UTC)) * R_x(-IERS::y_pole(Mjd_UTC));
    }
    
    
    //------------------------------------------------------------------------------
    //
    // Geodetic (class implementation)
    //
    // Purpose:
    //
    //   Class (with all public elements) for handling geodetic coordinates
    //
    //  p.185-189
    //  5.5 Geodetic Datums
    //------------------------------------------------------------------------------
    
    // Default constructor
    
    Geodetic::Geodetic ()
     : lon(0.0), lat(0.0), h(0.0)
    {
    }
    
    // Simple constructor
    
    Geodetic::Geodetic (double lambda, double phi, double alt)
    {
      lon=lambda; lat=phi; h=alt;
    }
    
    // Constructor for geodetic coordinates from given position
    
    Geodetic::Geodetic (Vector r,                         // Position vector [m]
                        double R_equ,                     // Equator radius [m]
                        double f)                         // Flattening
    {
    
      const double  eps     = 1.0e3*eps_mach;   // Convergence criterion 
      const double  epsRequ = eps*R_equ;
      const double  e2      = f*(2.0-f);        // Square of eccentricity
      
      const double  X = r(0);                   // Cartesian coordinates
      const double  Y = r(1);
      const double  Z = r(2);
      const double  rho2 = X*X + Y*Y;           // Square of distance from z-axis
      
      // Iteration 
    
      double  dZ, dZ_new, SinPhi;
      double  ZdZ, Nh, N;
    
      dZ = e2*Z;
      for(;;) {
        ZdZ    =  Z + dZ;
        Nh     =  sqrt ( rho2 + ZdZ*ZdZ ); 
        SinPhi =  ZdZ / Nh;                    // Sine of geodetic latitude
        N      =  R_equ / sqrt(1.0-e2*SinPhi*SinPhi); 
        dZ_new =  N*e2*SinPhi;
        if ( fabs(dZ-dZ_new) < epsRequ ) break;
        dZ = dZ_new;
      }
        
      // Longitude, latitude, altitude
    
      lon = atan2 ( Y, X );
      lat = atan2 ( ZdZ, sqrt(rho2) );
      h   = Nh - N;
    
    }
    
    // Position vector [m] from geodetic coordinates
    
    Vector Geodetic::Position (double R_equ,   // Equator radius [m]
                               double f    )   // Flattening
      const
    {  
    
      const double  e2     = f*(2.0-f);        // Square of eccentricity
      const double  CosLat = cos(lat);         // (Co)sine of geodetic latitude
      const double  SinLat = sin(lat);
    
      double  N;
      Vector  r(3);
          
      // Position vector 
    
      N = R_equ / sqrt(1.0-e2*SinLat*SinLat);
    
      r(0) =  (         N+h)*CosLat*cos(lon);
      r(1) =  (         N+h)*CosLat*sin(lon);
      r(2) =  ((1.0-e2)*N+h)*SinLat;
    
      return r;
    
    }
    
    // Transformation to local tangent coordinates
    
    Matrix Geodetic::LTC_Matrix () const
    {
      return LTCMatrix(lon,lat);
    }
    
    
    //------------------------------------------------------------------------------
    //
    // LTCMatrix
    //
    // Purpose:
    //
    //   Transformation from Greenwich meridian system to local tangent coordinates
    //
    // Input/Output:
    //
    //   lambda    Geodetic East longitude [rad]
    //   phi       Geodetic latitude [rad]
    //   <return>  Rotation matrix from the Earth equator and Greenwich meridian
    //             to the local tangent (East-North-Zenith) coordinate system
    //
    //------------------------------------------------------------------------------
    
    Matrix LTCMatrix (double lambda, double phi)
    {
      
      Matrix  M(3,3);
      double  Aux;
      
      // Transformation to Zenith-East-North System
      M = R_y(-phi)*R_z(lambda);
      
      // Cyclic shift of rows 0,1,2 to 1,2,0 to obtain East-North-Zenith system
      for (int j=0; j<3; j++) {
        Aux=M(0,j); M(0,j)=M(1,j); M(1,j)=M(2,j); M(2,j)= Aux;
      }
      
      return  M;
      
    }
    
    
    //------------------------------------------------------------------------------
    //
    // AzEl
    //
    // Purpose:
    //
    //   Computes azimuth and elevation from local tangent coordinates
    //
    // Input/Output:
    //
    //   s   Topocentric local tangent coordinates (East-North-Zenith frame)
    //   A   Azimuth [rad]
    //   E   Elevation [rad]
    //
    //------------------------------------------------------------------------------
    
    void AzEl ( const Vector& s, double& A, double& E ) 
    {
      A = atan2(s(0),s(1));
      A = ((A<0.0)? A+pi2 : A);
      E = atan ( s(2) / sqrt(s(0)*s(0)+s(1)*s(1)) );
    }
    
    
    //------------------------------------------------------------------------------
    //
    // AzEl
    //
    // Purpose:
    //
    //   Computes azimuth, elevation and partials from local tangent coordinates
    //
    // Input/Output:
    //
    //   s      Topocentric local tangent coordinates (East-North-Zenith frame)
    //   A      Azimuth [rad]
    //   E      Elevation [rad]
    //   dAds   Partials of azimuth w.r.t. s
    //   dEds   Partials of elevation w.r.t. s
    //
    //------------------------------------------------------------------------------
    
    void AzEl( const Vector& s, double& A, double& E, Vector& dAds, Vector& dEds ) 
    {
      const double rho = sqrt(s(0)*s(0)+s(1)*s(1));
      // Angles
      A = atan2(s(0),s(1));
      A = ((A<0.0)? A+pi2 : A);
      E = atan ( s(2) / rho );
      // Partials
      dAds = Vector( s(1)/(rho*rho), -s(0)/(rho*rho), 0.0 );
      dEds = Vector( -s(0)*s(2)/rho, -s(1)*s(2)/rho , rho ) / Dot(s,s);
    }
    
    
    Die Fehler-Ausgabe sieht so aus:
    Code:
    SAT_RefSys1.cpp: In function ‘Matrix EclMatrix(double)’:
    SAT_RefSys1.cpp:172: error: ‘R_x’ was not declared in this scope
    SAT_RefSys1.cpp: In function ‘Matrix PrecMatrix(double, double)’:
    SAT_RefSys1.cpp:217: error: ‘R_z’ was not declared in this scope
    SAT_RefSys1.cpp:217: error: ‘R_y’ was not declared in this scope
    SAT_RefSys1.cpp: In function ‘Matrix NutMatrix(double)’:
    SAT_RefSys1.cpp:437: error: ‘R_x’ was not declared in this scope
    SAT_RefSys1.cpp:437: error: ‘R_z’ was not declared in this scope
    SAT_RefSys1.cpp: In function ‘Matrix NutMatrixSimple(double)’:
    SAT_RefSys1.cpp:487: error: ‘R_x’ was not declared in this scope
    SAT_RefSys1.cpp:487: error: ‘R_z’ was not declared in this scope
    SAT_RefSys1.cpp: In function ‘Matrix GHAMatrix(double)’:
    SAT_RefSys1.cpp:620: error: ‘R_z’ was not declared in this scope
    SAT_RefSys1.cpp: In function ‘Matrix PoleMatrix(double)’:
    SAT_RefSys1.cpp:645: error: ‘R_y’ was not declared in this scope
    SAT_RefSys1.cpp:645: error: ‘R_x’ was not declared in this scope
    SAT_RefSys1.cpp: In function ‘Matrix LTCMatrix(double, double)’:
    SAT_RefSys1.cpp:773: error: ‘R_y’ was not declared in this scope
    SAT_RefSys1.cpp:773: error: ‘R_z’ was not declared in this scope
    make[4]: *** [SAT_RefSys1.o] Error 1
    
    Ich hoffe nicht, dass es zu weiteren Fehlern in den (anderen) Quellen kommt, da die urspr. Programmierer dieses Codes nicht mehr erreichbar sind.

    MfG,
    upperlimit
     
  6. rikola

    rikola Foren Gott

    Dabei seit:
    23.08.2005
    Beiträge:
    2.133
    Zustimmungen:
    0
    Wenn Du das Programm einfach nur kompilieren willst, ist das einfachste sicherlich, wenn Du Dir gcc-3.3.5 runterlaedst (gcc.gnu.org), kompilierst und damit dann das Programm kompilierst.

    Wenn Dir an dem Code was gelegen ist, musst Du die Fehlermeldungen Stueck fuer Stueck durcharbeiten und beheben. Suche mal nach der Funktion R_x und wo sie deklariert ist. Ich habe sie in den Dateien, die Du geschickt hast, nicht gefunden. Sie muss an der Stelle, an der sie in SAT_RefSys1.cpp zum ersten Mal benutzt wird, bereits bekannt, d.h., deklariert sein. In der SAT_VecMat.h ist sie als friend zur Klasse Matrix eingetragen, doch ist das keine Deklaration, soweit ich weiss.
     
  7. #6 Wolfgang, 29.06.2006
    Wolfgang

    Wolfgang Foren Gott

    Dabei seit:
    24.04.2005
    Beiträge:
    3.978
    Zustimmungen:
    0
    Ort:
    Erfurt
    Och ich kenne einige für gcc 3.3.X geschriebene Programme, die sich nicht so einfach mit gcc 4 compilieren lassen.

    Da sind eben Anpassungen nötig.
    Gruß Wolfgang
     
  8. rikola

    rikola Foren Gott

    Dabei seit:
    23.08.2005
    Beiträge:
    2.133
    Zustimmungen:
    0
    Schon klar. Ich bezog mich auf den von upperlimt angegebenen Code und v.a. die Fehlermeldung, die einfach typisch fuer eine fehlende Deklaration ist.
     
  9. Anzeige

    Vielleicht findest du HIER Antworten.
    Registrieren bzw. einloggen, um diese und auch andere Anzeigen zu deaktivieren
  10. #8 upperlimit, 29.06.2006
    upperlimit

    upperlimit Jungspund

    Dabei seit:
    11.04.2005
    Beiträge:
    16
    Zustimmungen:
    0
    Hallo,

    inzwischen habe ich den Code mit gcc 3.4.6 erfolgreich getestet, negativ mit gcc 4.1.1.


    Sofern rikola Recht hat, sollte es doch moeglich sein, die strengen ISO-Vorschriften vom gcc 4.1.0 mittels eines Schalters/Option zu entkraeften, so dass auch Code fuer aeltere gcc-Versionen problemlos kompiliert werden kann, richtig ?
    Ich habe bisher '-std=c++98', '-std=gnu++98' sowie '-ansi' ausprobiert, erfolglos.

    @rikola,
    Du hast recht, die Deklaration von R_x, R_y und R_z erfolgt in einer anderen Datei, genannt SAT_VecMat1.cpp (s. Code)

    Code:
    //------------------------------------------------------------------------------
    //
    // SAT_VecMat1.cpp
    // 
    // Purpose: 
    //
    //   Vector/matrix operations
    //   
    // Notes:
    //
    //   This software is protected by national and international copyright. 
    //   Any unauthorized use, reproduction or modificaton is unlawful and 
    //   will be prosecuted. Commercial and non-private application of the 
    //   software in any form is strictly prohibited unless otherwise granted
    //   by the authors.
    //
    //   The code is provided without any warranty; without even the implied 
    //   warranty of merchantibility or fitness for a particular purpose.
    //
    // Last modified:
    //
    //   2000/03/04  OMO  Final version (1st edition)
    //
    // (c) 1999-2000  O. Montenbruck, E. Gill
    //
    //------------------------------------------------------------------------------
    
    #include <cmath>
    #include <iostream>
    #include <iomanip>
    
    #include "SAT_VecMat.h"
    
    /*#ifdef __GNUC__  // GNU C++ adaptation
    #include "GNU_iomanip.h"
    #endif*/
    
    using namespace std;
    
    
    //------------------------------------------------------------------------------
    // 
    // LU_Decomp
    //
    // Purpose:
    //
    //   LU-Decomposition.
    //
    //   Given an nxn matrix A, this routine replaces it by the LU decomposition
    //   of a rowwise permutation of itself. A is output, arranged as in 
    //   equation (2.3.14) of Press et al. (1986); Indx is an ouput vector which
    //   records the row permutation effected by partial pivoting. This routine is 
    //   used in combination with LU_BackSub to solve linear equations or invert 
    //   a matrix.
    //
    // Input/output:
    //
    //   A       Square matrix; replaced by LU decomposition of permutation of A
    //           on output
    //   Indx    Permutation index vector
    //
    // Note:
    //
    //   Adapted from LUDCMP of Press et al. (1986).
    //
    //------------------------------------------------------------------------------
    
    void LU_Decomp ( Matrix& A, Vector& Indx )
    {
    
      // Constants
          
      const int    n    = A.size1();
      const double tiny = 1.0e-20;       // A small number
    
      // Variables
          
      int     i,j,imax,k;
      double  aAmax, Sum, Dum;
      Vector  V(n); 
    
      // Loop over rows to get scaling information
    
      for (i=0; i<n; i++) {
        aAmax = 0.0;
        for (j=0;j<n;j++) if (fabs(A(i,j)) > aAmax ) aAmax=fabs(A(i,j));
        if (aAmax==0.0) {
          // No nonzero largest element
          cerr << "ERROR: Singular matrix A in LU_Decomp";
          exit(1);
        };
        V(i) = 1.0/aAmax;           // V stores the implicit scaling of each row
      };
    
      // Loop over columns of Crout's method
          
      for ( j=0; j<n; j++ ) {
          
        if (j > 0) {
          for ( i=0; i<j; i++ ) {   // This is equation 2.3.12 except for i=j
            Sum = A(i,j);
            if (i>0) {
              for ( k=0; k<i; k++ )  Sum -= A(i,k)*A(k,j);
              A(i,j) = Sum;
            };
          };
        };
          
        aAmax=0.0;                  // Initialize for the search of the largest
                                    // pivot element
          
        for ( i=j; i<n; i++ ) {     // This is i=j of equation 2.3.12 and 
          Sum = A(i,j);             // i=j+1..N of equation 2.3.13
          if (j > 0) {
            for ( k=0; k<j; k++ ) Sum -= A(i,k)*A(k,j);
            A(i,j) = Sum;
          };
          Dum = V(i)*fabs(Sum);     // Figure of merit for the pivot
          if (Dum >= aAmax) {       // Is it better than the best so far ?
            imax  = i;
            aAmax = Dum;
          };
        };
          
        if (j != imax) {            // Do we need to interchange rows?
          for ( k=0; k<n; k++) {    // Yes, do so ...
            Dum = A(imax,k);
            A(imax,k) = A(j,k);
            A(j,k) = Dum;
          }
          V(imax) = V(j);           // Also interchange the scale factor 
        };
          
        Indx(j) = imax;
          
        if (j != n-1) {             // Now finally devide by the pivot element
          if (A(j,j) == 0.0) {      // If the pivot element is zero the matrix 
            A(j,j) = tiny;          // is singular (at least to the precision of
          };                        // the algorithm). For some applications on
          Dum=1.0/A(j,j);           // singular matrices, it is desirable to 
          for (i=j+1;i<n;i++) {     // substitude tiny for zero. 
            A(i,j)=A(i,j)*Dum;
          };
        };
    
      };   // Go back for the next column in the reduction
    
      if (A(n-1,n-1)==0.0) A(n-1,n-1)=tiny; 
    
    };
    
    
    //------------------------------------------------------------------------------
    //
    // LU_BackSub
    //
    // Purpose:
    //
    //   LU Backsubstitution
    //
    //   Solves the set of n linear equations Ax=b. Here A is input, not as the 
    //   matrix A but rather as its LU decomposition, determined by the function
    //   LU_Decomp. b is input as the right-hand side vector b, and returns with
    //   the solution vector x. A and Indx are not modified by this function and 
    //   can be left in place for successive calls with different right-hand 
    //   sides b. This routine takes into account the posssibility that B will  
    //   begin with many zero elements, so it is efficient for use in matrix
    //   inversions.
    //
    // Input/output:
    //
    //   A       LU decomposition of permutation of A
    //   Indx    Permutation index vector
    //   b       Right-hand side vector b; replaced by solution x of Ax=b on output
    //
    //------------------------------------------------------------------------------
    
    void LU_BackSub ( Matrix& A, Vector& Indx, Vector& b )
    { 
    
      // Constants
          
      const int  n = A.size1();
    
      // Local variables
    
      int     ii,i,ll,j;
      double  Sum;
    
      //
      // Start
      //
    
      ii = -1;                      // When ii is set to a nonegative value, it will
                                    // become the first nonvanishing element of B. 
      for (i=0; i<n; i++) {         // We now do the forward substitution.
        ll = (int) Indx(i);         // The only wrinkle is to unscramble the 
        Sum = b(ll);                // permutation as we go.
        b(ll) = b(i);
        if (ii != -1) {
          for (j=ii; j<i; j++) Sum -= A(i,j)*b(j);
        }
        else {
          if (Sum != 0.0) ii = i;   // A nonzero element was encountered, so from 
        };                          // now on we will have to do the sums in the
        b(i) = Sum;                 // loop above.
       };
          
       for (i=n-1; i>=0; i--) {     // Now we do the backsubstitution, eqn 2.3.7.
         Sum=b(i);
         if (i<n-1) {
           for (j=i+1;j<n;j++) {
              Sum = Sum-A(i,j)*b(j);
           };
         };
         b(i) = Sum/A(i,i);         // Store a component of the solution vector X.
       };
    
    };
    
    
    
    //------------------------------------------------------------------------------
    //
    // Vector (class implementation)
    //
    // Purpose:
    //
    //   Vector data type and associated operations
    //
    //------------------------------------------------------------------------------
    
    
    // Constructors, destructor
    
    Vector::Vector ()                              // Vector without elements
      : n(0) 
    {
      v = 0;
    }
    
    Vector::Vector (int Size)                      // Creates null vector
     : n(Size)
    {
      v = new double [Size];
      for (int i=0; i<Size; i++) v[i]=0.0;
    }
    
    Vector::Vector (const Vector& V)               // Vector copy
     : n(V.n)
    {
      v = new double [V.n];
      for (int i=0; i<V.n; i++) v[i]=V.v[i];
    }
    
    Vector::Vector (const double* p, int N)        // Array copy
      : n(N)
    {
      v = new double [N];
      for (int i=0; i<N; i++) v[i]=p[i];
    }
    
    Vector::Vector (double x, double y, double z)  // 3dim-Vector
      : n(3)
    {
      v = new double [3];
      v[0]=x; v[1]=y; v[2]=z;
    }
    
    Vector::Vector (double x, double y, double z,   // 6dim-Vector
                    double X, double Y, double Z)  
      : n(6)
    {
      v = new double [6];
      v[0]=x; v[1]=y; v[2]=z;
      v[3]=X; v[4]=Y; v[5]=Z;
    }
    
    Vector::~Vector() 
    { 
      delete [] v; 
    }
    
    // Component access
    
    Vector Vector::slice (int first, int last) const
    {
      Vector Aux(last-first+1);
      for (int i=first; i<=last; i++) Aux.v[i-first]=v[i];
      return Aux;
    }
    
    
    // Square root of vector elements
    
    Vector Vector::Sqrt()
    {
      Vector Aux(n);
      for (int i=0; i<n; i++) Aux.v[i]=sqrt(v[i]);
      return Aux;
    }
    
    
    // Assignment
    
    Vector& Vector::operator=(const double value)
    {
      for (int i=0; i<n; i++) v[i]=value;
      return (*this);
    }
    
    Vector& Vector::operator=(const Vector& V)
    {
      if (this == &V) return (*this);
      // Allocate vector if still empty
      if (v==0) {
        n = V.n; 
        v = new double [V.n];
      };
      // Check dimension
      if (n!=V.n) {
        cerr << "ERROR: Incompatible sizes in Vector operator=(Vector)" << endl;
        exit(1);
      };
      // Copy elements
      for (int i=0; i<n; i++) v[i]=V.v[i];
      return (*this);
    }
    
    
    // Concatenation
    
    Vector Stack (const Vector& a, const Vector& b)
    {
      int    i;
      Vector c(a.size()+b.size());
      for (i=0;i<a.size();i++) c(i)=a(i);
      for (i=0;i<b.size();i++) c(i+a.size())=b(i);
      return c;
    }
    
    
    // Vector from polar angles
    
    Vector VecPolar (double azim, double elev, double r)
    {
      return Vector(r*cos(azim)*cos(elev),r*sin(azim)*cos(elev),r*sin(elev));
    }
    
    
    // Vector addition/subtraction with assignment
    
    void Vector::operator += (const Vector& V)
    {
      if (n!=V.n) {
        cerr << "ERROR: Incompatible shape in Vector operator+=(Vector)" << endl;
        exit(1);
      };
      for (int i=0; i<n; i++) v[i]+=V.v[i];
    }
    
    void Vector::operator -= (const Vector& V)
    {
      if (n!=V.n) {
        cerr << "ERROR: Incompatible shape in Vector operator-=(Vector)" << endl;
        exit(1);
      };
      for (int i=0; i<n; i++) v[i]-=V.v[i];
    }
    
    
    // Dot product, norm, cross product
    
    double Dot (const Vector& left, const Vector& right)
    {
      if (left.n!=right.n) {
        cerr << "ERROR: Incompatible shape in Dot(Vector,Vector)" << endl;
        exit(1);
      };
      double Sum = 0.0;
      for (int i=0; i<left.n; i++) Sum+=left.v[i]*right.v[i];
      return Sum;
    }
    
    double Norm (const Vector& V)
    {
      return sqrt(Dot(V,V));
    }
    
    Vector Cross (const Vector& left, const Vector& right)
    {
      if ( (left.n!=3) || (right.n!=3) ) {
        cerr << "ERROR: Invalid dimension in Cross(Vector,Vector)" << endl;
        exit(1);
      };
      Vector Result(3);
      Result.v[0] = left.v[1]*right.v[2] - left.v[2]*right.v[1];
      Result.v[1] = left.v[2]*right.v[0] - left.v[0]*right.v[2];
      Result.v[2] = left.v[0]*right.v[1] - left.v[1]*right.v[0];
      return Result;
    }
    
    
    // Scalar multiplication and division of a vector
    
    Vector operator * (double value, const Vector& V)
    {
      Vector Aux(V.n);
      for (int i=0; i<V.n; i++) Aux.v[i]=value*V.v[i];
      return Aux;
    }
    
    Vector operator * (const Vector& V, double value)
    {
      return value*V;
    }
    
    Vector operator / (const Vector& V, double value)
    {
      Vector Aux(V.n);
      for (int i=0; i<V.n; i++) Aux.v[i]=V.v[i]/value;
      return Aux;
    }
    
    
    // Negation of a vector (unary minus)
    
    Vector operator - (const Vector& V)
    {
      Vector Aux(V.n);
      for (int i=0; i<V.n; i++) Aux.v[i]=-V.v[i];
      return Aux;
    }
    
    
    // Vector addition and subtraction
    
    Vector operator + (const Vector& left, const Vector& right)
    {
      if (left.n!=right.n) {
        cerr << "ERROR: Incompatible shape in +(Vector,Vector)" << endl;
        exit(1);
      };
      Vector Aux(left.n);
      for (int i=0; i<left.n; i++) Aux.v[i]=left.v[i]+right.v[i];
      return Aux;
    }
    
    Vector operator - (const Vector& left, const Vector& right)
    {
      if (left.n!=right.n) {
        cerr << "ERROR: Incompatible shape in -(Vector,Vector)" << endl;
        exit(1);
      };
      Vector Aux(left.n);
      for (int i=0; i<left.n; i++) Aux.v[i]=left.v[i]-right.v[i];
      return Aux;
    }
    
    // Vector output
    
    ostream& operator << (ostream& os, const Vector& Vec)
    {
      int w = os.width();
      for (int i=0; i<Vec.size(); i++)
        os << setw(w) << Vec(i);
    //os << endl;
      return os;
    }
    
    
    
    //------------------------------------------------------------------------------
    //
    // Matrix (class implementation)
    //
    // Purpose:
    //
    //   Matrix data type and associated operations
    //
    //------------------------------------------------------------------------------
    
    
    // Constructors
    
    Matrix::Matrix ()                          // Matrix without elements
      : n(0), m(0) 
    {
      M = 0;
    }
    
    Matrix::Matrix (int dim1, int dim2)        // Nullmatrix of specified shape 
      : n(dim1), m(dim2)
    {
      int i,j;
      // Memory allocation
      M = new double*[dim1];
      for (i=0; i<dim1; i++) M[i] = new double[dim2];
      // Initialization
      for (i=0; i<dim1; i++) {
        for (j=0; j<dim2; j++) M[i][j]=0.0;
      }
    }
    
    Matrix::Matrix (const Matrix& M_)          // Copy
    {
      int i,j;
      n = M_.n;
      m = M_.m;
      // Memory allocation
      M = new double*[n];
      for (i=0; i<n; i++) M[i] = new double[m];
      // Initialization
      for (i=0; i<n; i++) {
        for (j=0; j<m; j++) M[i][j]=M_.M[i][j];
      }
    }
    
    Matrix::Matrix (const double* p, int dim1, int dim2)   // Array copy
    {
      int i,j;
      n = dim1;
      m = dim2;
      // Memory allocation
      M = new double*[n];
      for (i=0; i<n; i++) M[i] = new double[m];
      // Initialization
      for (i=0; i<n; i++) {
        for (j=0; j<m; j++) M[i][j]=p[i*dim2+j];
      }
      
    }
    
    Matrix::~Matrix() 
    { 
      for (int i=0; i<n; i++) delete[] M[i];
      delete [] M; 
    };
    
    
    // Assignment
    
    Matrix& Matrix::operator=(const double value)
    {
      for (int i=0; i<n; i++) 
        for (int j=0; j<n; j++) 
          M[i][j]=value;
      return (*this);
    }
    
    Matrix& Matrix::operator=(const Matrix& M_)
    {
      int i,j;
      if (this == &M_) return (*this);
      // Allocate matrix if still empty
      if (M==0) {
        n = M_.n; 
        m = M_.m;
        M = new double*[n];
        for (i=0; i<n; i++) M[i] = new double[m];
      };
      if ( (n!=M_.n) || (m!=M_.m) ) {
        cerr << "ERROR: Incompatible shapes in Matrix operator=(Matrix)" << endl;
        exit(1);
      };
      for (i=0; i<n; i++) 
        for (j=0; j<m; j++) 
          M[i][j]=M_.M[i][j];
      return (*this);
    }
    
    
    // Component access
    
    Vector Matrix::Col(int j) const
    {
      Vector Res(n);
      for (int i=0; i<n; i++)  Res.v[i]=M[i][j];
      return Res;
    }
    
    Vector Matrix::Row(int i) const
    {
      Vector Res(m);
      for (int j=0; j<m; j++)  Res.v[j]=M[i][j];
      return Res;
    }
    
    Vector Matrix::Diag() const
    {
      if (n!=m) {
        cerr << "ERROR: Invalid shape in Matrix.Diag()" << endl;
        exit(1);
      };
      Vector Vec(n);
      for (int i=0; i<n; i++) Vec.v[i] = M[i][i];
      return Vec;
    }
    
    
    Matrix Matrix::slice(int first_row, int last_row, int first_col, int last_col)
    {
      if (first_row<0 || last_row<first_row || n-1<last_row ||
          first_col<0 || last_col<first_col || m-1<last_col) {
        cerr << "ERROR: Invalid arguments in Matrix.slice()" << endl;
        exit(1);
      };
      Matrix Aux(last_row-first_row+1,last_col-first_col+1);
      for (int i=0;i<=last_row-first_row;i++)
        for (int j=0;j<=last_col-first_col;j++)
           Aux(i,j) = M[i+first_row][j+first_col];
      return Aux;
    }
    
    
    void Matrix::SetCol(int j, const Vector& Col) 
    {
      if (Col.size()!=n) {
        cerr << "ERROR: Incompatible shapes in Matrix.SetCol()" << endl;
        exit(1);
      };
      if (j<0 || m<=j) {
        cerr << "ERROR: Column index out of range in Matrix.SetCol()" << endl;
        exit(1);
      };
      for (int i=0; i<n; i++) M[i][j]=Col(i);
    }
    
    void Matrix::SetRow(int i, const Vector& Row)
    {
      if (Row.size()!=m) {
        cerr << "ERROR: Incompatible shapes in Matrix.SetRow()" << endl;
        exit(1);
      };
      if (i<0 || n<=i) {
        cerr << "ERROR: Row index out of range in Matrix.SetRow()" << endl;
        exit(1);
      };
      for (int j=0; j<n; j++) M[i][j]=Row(j);
    }
    
    
    // Unit matrix
    
    Matrix Id(int Size)
    {
      Matrix Aux(Size,Size);      
      for (int i=0; i<Size; i++) Aux.M[i][i] = 1.0;
      return Aux;
    }
    
    
    // Diagonal matrix
    
    Matrix Diag(const Vector& Vec)
    {
      Matrix Mat(Vec.n,Vec.n);
      for (int i=0; i<Vec.n; i++) Mat.M[i][i] = Vec.v[i];
      return Mat;
    }
    
    
    // Elementary rotations
    
    Matrix R_x(double Angle)
    {
      const double C = cos(Angle);
      const double S = sin(Angle);
      Matrix U(3,3);
      U.M[0][0] = 1.0;  U.M[0][1] = 0.0;  U.M[0][2] = 0.0;
      U.M[1][0] = 0.0;  U.M[1][1] =  +C;  U.M[1][2] =  +S;
      U.M[2][0] = 0.0;  U.M[2][1] =  -S;  U.M[2][2] =  +C;
      return U;
    }
    
    Matrix R_y(double Angle)
    {
      const double C = cos(Angle);
      const double S = sin(Angle);
      Matrix U(3,3);
      U.M[0][0] =  +C;  U.M[0][1] = 0.0;  U.M[0][2] =  -S;
      U.M[1][0] = 0.0;  U.M[1][1] = 1.0;  U.M[1][2] = 0.0;
      U.M[2][0] =  +S;  U.M[2][1] = 0.0;  U.M[2][2] =  +C;
      return U;
    }
    
    Matrix R_z(double Angle)
    {
      const double C = cos(Angle);
      const double S = sin(Angle);
      Matrix U(3,3);
      U.M[0][0] =  +C;  U.M[0][1] =  +S;  U.M[0][2] = 0.0;
      U.M[1][0] =  -S;  U.M[1][1] =  +C;  U.M[1][2] = 0.0;
      U.M[2][0] = 0.0;  U.M[2][1] = 0.0;  U.M[2][2] = 1.0;
      return U;
    }
    
    
    // Transposition
    
    Matrix Transp(const Matrix& Mat)
    {
      Matrix T(Mat.m,Mat.n);
      for ( int i=0; i<T.n; i++ )
        for ( int j=0; j<T.m; j++ )
          T.M[i][j] = Mat.M[j][i];
      return T;
    }
    
    
    // Inverse
    
    Matrix Inv(const Matrix& Mat)
    {
      const int n = Mat.n;
    
      Matrix LU(n,n), Inverse(n,n);
      Vector b(n), Indx(n);
    
      if (Mat.m!=Mat.n) {
        cerr << "ERROR: Invalid shape in Inv(Matrix)" << endl;
        exit(1);
      };
    
      // LU decomposition 
    
      LU = Mat;
      LU_Decomp ( LU, Indx );
    
      // Solve Ax=b for  unit vectors b_1..b_n
    
      for (int j=0; j<n; j++ ) {
        b=0.0; b(j)= 1.0;                     // Set b to j-th unit vector 
        LU_BackSub ( LU, Indx, b );           // Solve Ax=b 
        Inverse.SetCol(j,b);                  // Copy result
      };
    
      return Inverse;
    
    }
    
    
    // Scalar multiplication and division of a vector
    
    Matrix operator * (double value, const Matrix& Mat)
    {
      Matrix Aux(Mat.n,Mat.m);
      for (int i=0; i<Mat.n; i++) 
        for (int j=0; j<Mat.m; j++) 
          Aux.M[i][j]=value*Mat.M[i][j];
      return Aux;
    }
    
    Matrix operator * (const Matrix& Mat, double value)
    {
      return value*Mat;
    }
    
    Matrix operator / (const Matrix& Mat, double value)
    {
      Matrix Aux(Mat.n,Mat.m);
      for (int i=0; i<Mat.n; i++) 
        for (int j=0; j<Mat.m; j++) 
          Aux.M[i][j]=Mat.M[i][j]/value;
      return Aux;
    }
    
    
    // Unary minus
    
    Matrix operator - (const Matrix& Mat)
    {
      Matrix Aux(Mat.n,Mat.m);
      for (int i=0; i<Mat.n; i++) 
        for (int j=0; j<Mat.m; j++) 
          Aux.M[i][j]=-Mat.M[i][j];
      return Aux;
    }
    
    
    // Matrix addition and subtraction
    
    Matrix operator + (const Matrix& left, const Matrix& right)
    {
      if ( (left.n!=right.n) || (left.m!=right.m) ) {
        cerr << "ERROR: Incompatible shape in +(Matrix,Matrix)" << endl;
        exit(1);
      };
      Matrix Aux(left.n,left.m);
      for (int i=0; i<left.n; i++) 
        for (int j=0; j<left.m; j++) 
          Aux.M[i][j] = left.M[i][j] + right.M[i][j];
      return Aux;
    }
    
    Matrix operator - (const Matrix& left, const Matrix& right)    
    {
      if ( (left.n!=right.n) || (left.m!=right.m) ) {
        cerr << "ERROR: Incompatible shape in -(Matrix,Matrix)" << endl;
        exit(1);
      };
      Matrix Aux(left.n,left.m);
      for (int i=0; i<left.n; i++) 
        for (int j=0; j<left.m; j++) 
          Aux.M[i][j] = left.M[i][j] - right.M[i][j];
      return Aux;
    }
    
    
    // Matrix product
    
    Matrix operator * (const Matrix& left, const Matrix& right)
    {
      if (left.m!=right.n) {
        cerr << "ERROR: Incompatible shape in *(Matrix,Matrix)" << endl;
        exit(1);
      };
      Matrix Aux(left.n,right.m);
      double Sum;
      for (int i=0; i<left.n; i++) 
        for (int j=0; j<right.m; j++) {
          Sum = 0.0;
          for (int k=0; k<left.m; k++) 
            Sum += left.M[i][k] * right.M[k][j];
          Aux.M[i][j] = Sum;
        }
      return Aux;
    }
    
    
    // Vector/matrix product
    
    Vector operator * (const Matrix& Mat, const Vector& Vec)
    {
      if (Mat.m!=Vec.n) {
        cerr << "ERROR: Incompatible shape in *(Matrix,Vector)" << endl;
        exit(1);
      };
      Vector Aux(Mat.n);
      double Sum;
      for (int i=0; i<Mat.n; i++) {
        Sum = 0.0;
        for (int j=0; j<Mat.m; j++) 
          Sum += Mat.M[i][j] * Vec.v[j];
        Aux.v[i] = Sum;
      }
      return Aux;
    }
    
    Vector operator * (const Vector& Vec, const Matrix& Mat)
    {
      if (Mat.n!=Vec.n) {
        cerr << "ERROR: Incompatible shape in *(Vector,Matrix)" << endl;
        exit(1);
      };
      Vector Aux(Mat.m);
      double Sum;
      for (int j=0; j<Mat.m; j++) {
        Sum = 0.0;
        for (int i=0; i<Mat.n; i++) 
          Sum += Vec.v[i] * Mat.M[i][j];
        Aux.v[j] = Sum;
      }
      return Aux;
    }
    
    
    // Dyadic product
    
    Matrix Dyadic (const Vector& left, const Vector& right)
    {
      Matrix Mat(left.n,right.n);
      for (int i=0;i<left.n;i++)
        for (int j=0;j<right.n;j++)
          Mat.M[i][j] = left.v[i]*right.v[j];
      return Mat;
    }
    
    
    // Matrix output
    
    ostream& operator << (ostream& os, const Matrix& Mat)
    {
      int w = os.width();
      for (int i=0; i<Mat.size1(); i++) {
        for (int j=0; j<Mat.size2(); j++)
          os << setw(w) << Mat(i,j);
        os << endl;
      }
      return os;
    }
    
    

    MfG,
    upperlimit
     
  11. #9 mindless, 29.06.2006
    mindless

    mindless Grünschnabel

    Dabei seit:
    29.06.2006
    Beiträge:
    3
    Zustimmungen:
    0
    Ort:
    Hamburg
    naja mit der option -std=c++98 kannst du ja auch nicht viel erreichen.
    Damit sagst du ja du willst den c++ standard von 98, dies ist der aktuelle
    ISO-standard. Den willst du ja grade nicht!
    um ehrlich zu sein, kann ich dir leider auch keinen schalter nennen, der den
    standard "abschlatet", btw denke ich auch, dass dies vll nicht wirklich der
    richtige weg ist.
     
Thema:

gcc Versionsunterschied -> Problem

Die Seite wird geladen...

gcc Versionsunterschied -> Problem - Ähnliche Themen

  1. Virtualboxproblem

    Virtualboxproblem: Hallo Folgnedes Problem: Ich kann ohne Probleme In Debian-Sid, 64Bit und virtualbox-dkms Unixoide installieren, wenn es sich um eine...
  2. grub-pc Probleme bei upgrade

    grub-pc Probleme bei upgrade: Hallo, ich habe beim dist-upgrade folgendes Problem: ---------- Nach dieser Operation werden 0 B Plattenplatz zusätzlich benutzt. Trigger für...
  3. Problem bei apt-get upgrade (Kali 2.0)

    Problem bei apt-get upgrade (Kali 2.0): Hi, seit paar Tagen habe ich einige Probs. Dachte es sei mal an der Zeit für ein Update und nun tauchen da einige Fehler auf: Es müssen noch 0 B...
  4. Probleme mit YUM

    Probleme mit YUM: Hallo, ich habe CentOs 7 als Dualboot mit Windows 7 auf einen Dell Latitude E5510 installiert. Dies hat soweit auch alles geklappt. Leider habe...
  5. Problem mit Win-Zugriff auf SAMBA

    Problem mit Win-Zugriff auf SAMBA: Hallo rundherum, vielleicht hat jemand eine Idee... Habe jetzt schon 1 1/2 Tage damit verbraten und den Fehler nicht gefunden. Problem:...