试图模拟双摆的分段错误

Question

我正在尝试编写一个程序来模拟双摆，很像 this program。当我运行此代码出现段错误时，但我找不到导致它的原因。


#include <iostream>
#include <cmath>
#include <vector>

#define PI 3.14159265 // Pi
#define G 9.8 // Acceleration due to gravity
#define L1 1.0 // Pendulum 1 Length ( meters )
#define M1 1.0 // Pendulum 1 Mass ( kilograms )
#define L2 1.0 // Pendulum 2 Length ( meters )
#define M2 1.0 // Pendulum 2 Mass ( kilograms )
#define N 4 // Number of equations. Don't change

using namespace std;

void calculateDerivatives ( vector<float> &dydx , vector<float> &yIn )
{
   // Calculate derivatives and put them in dydx vector
   dydx.clear();
   float temp1 , temp2 , delta;

   delta = yIn[2] - yIn[0];

   // Calculate sins and cosines to lower computation cost
   float sD = sin ( delta );
   float cD = cos ( delta );
   float s2D = sin ( 2 * delta );

   dydx.push_back ( yIn[1] );

   temp1 = ( M2 * L1 * pow ( yIn[1] , 2 ) * (1 / 2) * s2D );
   temp1 += ( M2 * G * sin ( yIn[2] ) * cD );
   temp1 += ( M2 * L2 * pow ( yIn[3] , 2 ) * sD );
   temp1 -= ( ( M1 + M2 ) * G * sin ( yIn[0] ) ) / ( ( M1 + M2 ) * L1 - M2 *L1 * pow ( cD , 2 ) );
   dydx.push_back ( temp1 );

   dydx.push_back ( yIn[3] );

   temp2 = ( -M2 * L2 * pow ( yIn[3] , 2 ) * ( 1 / 2 ) * s2D );
   temp2 += ( M1 + M2 ) * G * sin ( yIn[0] ) * cD;
   temp2 -= ( M1 + M2 ) * L1 * pow ( yIn[1] , 2) * sD;
   temp2 -= ( ( M1 + M2 ) * G * sin ( yIn[2] ) ) / ( ( L2 / L1 ) * delta );
   dydx.push_back ( temp2 );
}

void doRungeKutta ( const float x , const float h , vector<float> &yIn , vector<float> &yOut )
{
   // Use the Runge Kutta method to solve the DE
   vector<float> dydx;
   vector<float> dydxt;
   vector<float> yt;
   vector< vector<float> > kuttas;

   // Step 1
   calculateDerivatives ( dydx , yIn );
   vector<float> temp;
   for ( int i = 0; i < N; i++ )
   {
      temp.push_back( h * dydx[i] );
      yt.push_back ( yIn[i] + ( 0.5 ) * temp[i] );
   }
   kuttas.push_back ( temp );

   // Step 2
   calculateDerivatives ( dydxt , yt );
   for ( int i = 0; i < N; i++ )
   {
      temp[i] = ( h * dydxt[i] );
      yt[i] = ( yIn[i] + ( 0.5 ) * temp[i] );
   }
   kuttas.push_back ( temp );

   // Step 3
   calculateDerivatives ( dydxt , yt );
   for ( int i = 0; i < N; i++ )
   {
      temp[i] = ( h * dydxt[i] );
      yt[i] = ( yIn[i] + temp[i] );
   }
   kuttas.push_back ( temp );

   // Step 4
   calculateDerivatives ( dydxt , yt );
   for ( int i = 0; i < N; i++ )
   {
      temp[i] = ( h * dydxt[i] );
      yOut.push_back ( yIn[i] + kuttas[0][i] / 6.0 + kuttas[1][i] / 3.0 + kuttas[2][i] / 3.0 + kuttas[3][i] / 6.0 );
   }

}

int main( int argc , char *argv[] )
{
   int steps;
   float time_min , time_max , sTheta1 , sTheta2 , sAngVel1 , sAngVel2;
   float degRad = PI / 180;
   vector<float> yIn , yOut , times , theta1 , theta2 , angVel1 , angVel2;

   // Get command line inputs
   time_min = atof ( argv[1] );
   time_max = atof ( argv[2] );
   sTheta1 = atof ( argv[3] );
   sAngVel1 = atof ( argv[4] );
   sTheta2 = atof ( argv[5] );
   sAngVel2 = atof ( argv[6] );
   steps = atoi ( argv[7] );

   int h = ( time_max - time_min ) / ( steps - 1 );
   for ( int i = 0; i < steps; i++ )
   {
      times[i] = time_min + h * i;
   }

   theta1.push_back( sTheta1 * degRad );
   theta2.push_back( sTheta2 * degRad );
   angVel1.push_back( sAngVel1 );
   angVel2.push_back( sAngVel2 );

   cout << times[0] << " " << theta1[0] << " "
        << angVel1[0] << " " << theta2[0] << " " << angVel2[0];

   for ( int i = 0; i < ( steps - 1 ); i++ )
   {
      yIn.push_back( theta1[i] );
      yIn.push_back( angVel1[i] );
      yIn.push_back( theta2[i] );
      yIn.push_back( angVel2[i] );

      doRungeKutta ( times[i] , h , yIn , yOut );

      theta1.push_back( yOut[0] );
      angVel1.push_back( yOut[1] );
      theta2.push_back( yOut[2] );
      angVel2.push_back( yOut[3] );
      cout << times[i + 1] << " " << theta1[i + 1] << " "
           << angVel1[i + 1] << " " << theta2[i + 1] << " " << angVel2[i + 1];

   }

   return 0;
}

抱歉，太长了，我只是不知道是什么原因造成的，所以我无法将其缩短为有问题的部分。我一遍又一遍地检查代码，寻找我在声明数组或越界访问时可能犯的某种错误，但我找不到我的错误。

运行 gdb 我明白了：

$ g++ dp.cpp
$ ./a.out 0 10 90 0 90 0 1000
Segmentation fault (core dumped)
$ gdb a.out
GNU gdb (Ubuntu 8.1-0ubuntu3.2) 8.1.0.20180409-git
Copyright (C) 2018 Free Software Foundation, Inc.
License GPLv3+: GNU GPL version 3 or later <http://gnu.org/licenses/gpl.html>
This is free software: you are free to change and redistribute it.
There is NO WARRANTY, to the extent permitted by law.  Type "show copying"
and "show warranty" for details.
This GDB was configured as "x86_64-linux-gnu".
Type "show configuration" for configuration details.
For bug reporting instructions, please see:
<http://www.gnu.org/software/gdb/bugs/>.
Find the GDB manual and other documentation resources online at:
<http://www.gnu.org/software/gdb/documentation/>.
For help, type "help".
Type "apropos word" to search for commands related to "word"...
Reading symbols from a.out...(no debugging symbols found)...done.
(gdb) run 0 10 90 0 90 0 1000
Starting program: a.out 0 10 90 0 90 0 1000

Program received signal SIGSEGV, Segmentation fault.
0x0000000008001743 in doRungeKutta(float, float, std::vector<float, std::allocator<float> >&, std::vector<float, std::allocator<float> >&) ()
(gdb) where
#0  0x0000000008001743 in doRungeKutta(float, float, std::vector<float, std::allocator<float> >&, std::vector<float, std::allocator<float> >&) ()
#1  0x0000000008001c8a in main ()
(gdb)

这对我没有帮助。我尝试在doRungeKutta函数中搜索找到问题所在，可惜我找不到。

编辑：使用向量而不是分配内存。但是，相同的 gdb 输出仍然会出现段错误。

Answer 1

这里有一个错误：

void doRungeKutta ( const float x , const float h , vector<float> &yIn , vector<float> &yOut )
{
   //...
   vector< vector<float> > kuttas;

   // Step 1
   calculateDerivatives ( dydx , yIn );
   //...
   kuttas.push_back ( temp );  // Here is the first one

   // Step 2
   calculateDerivatives ( dydxt , yt );
   //.. 
   kuttas.push_back ( temp );  // Here is the second one

   // Step 3
   calculateDerivatives ( dydxt , yt );
   //...   
   kuttas.push_back ( temp ); // Here is the third one

   // Step 4
   calculateDerivatives ( dydxt , yt );
   for ( int i = 0; i < N; i++ )
   {
      temp[i] = ( h * dydxt[i] );
      yOut.push_back ( yIn[i] + kuttas[0][i] / 6.0 + kuttas[1][i] / 3.0 + kuttas[2][i] 
                  / 3.0 + kuttas[3][i] / 6.0 );  // <-- Error -- there is no kutta[3]
   }
}

kuttas 向量的大小为 3，但您正在尝试访问 kutta[3]，这是越界访问。因此内存覆盖和未定义的行为。

您还有其他错误，例如某些 vector 项具有 -inf 值，但明显的错误是指出的错误。

试图模拟双摆的分段错误

Segmentation Fault Trying To Simulate Double Pendulum

c++

gdb

segmentation-fault