c - C语言中ceil()和floor()的输出为奇数

Question

我正在做功课来实现一个 C 程序，在我的工作过程中，有一个步骤是我需要获取 double 型数字的整数部分。所以我选择ceil()或floor()来实现这一点。但输出是不可预测的，并且与预期相去甚远。

整个程序如下:

 /*
 ************************************************************************
  Includes
 ************************************************************************
*/
#include <stdio.h>
#include <stdlib.h>
#include <gsl/gsl_rng.h>
#include <time.h>
#include <math.h>

/* Solvent Particle Propersities*/
typedef struct 
  {
    double vx,vy,rx,ry ;  /*  velocity AND position               */
    double    cell;                /*  index of grid  */
   }solvent_p;

  /* Cell Propersities*/
 typedef struct 
  {
  double vcm_x,vcm_y ;  /*  center of mass velocity  */
  int    number;        /*  number of solvent in the cell  */
  double roation_angle; /*  roation_angle of the cell */
  }cell_p;

/* periodic boundary condition judging and adjusting fuction PBC */

 void PBC(solvent_p *sol)
 {
  double L = 20.0;  // box size 20

  if(sol->rx >20) sol->rx=sol->rx-L;
  if(sol->rx < 0) sol->rx=sol->rx+L;
  if(sol->ry >20) sol->ry=sol->ry-L;
  if(sol->ry < 0) sol->ry=sol->ry+L;

 }

int main(int argc, char **argv)
{

   // Randome setup generates random numbers from GSL functions 
   const gsl_rng_type * T;
   gsl_rng * r;
   T = gsl_rng_default;
   gsl_rng_default_seed = ((unsigned long)(time(NULL))); //设seed值为当前时间 
   r = gsl_rng_alloc (T);

   solvent_p solvent[4000];
   int i,j,k,index=0;
   cell_p grid[400];
   double alpha=90.0; //roation angle

  /*  Iniinitializing solvent  
   *  All vx=0
   *  half vy = sqrt(3) and half vy=-sqrt(3) to make momentum zero and another requirement is the overall energy is 6000 equals energy of temperature=1 with 4000 solvent 3NkT/2 ,assuming mass of solvent = 1 ,which is all a test quantity 
  *  rx and ry are random number generated by GSL library
  *  cell=20*(ry+rx) is an index of which cell the solvent is in
  */
  for(i=0;i<=3999;i++)
  {
         solvent[i].vx=0.0;

         if(i<=1999)
              solvent[i].vy=sqrt(3);
         else
              solvent[i].vy=-sqrt(3);


      solvent[i].rx =20.0 * gsl_rng_uniform_pos(r);

      solvent[i].ry =20.0 * gsl_rng_uniform_pos(r);

     //printf("%f \t %f \n",solvent[i].rx,solvent[i].ry);

    solvent[i].cell=20*(floor(solvent[i].ry))+floor(solvent[i].rx)+1;
    }

     // grid setting up
     for (i=0;i<=399;i++)
   {
        grid[i].vcm_x=0;
        grid[i].vcm_y=0;
        grid[i].number=0;
        grid[i].roation_angle=0.0;
   }


        /*  Begining Simulation Work
         * 
         *  Fist process is preparing the system to equilibrium for 10000 processes
          * 
          *  the whole process involving two steps steaming and collision and the two steps are conducted one by one in our simulation 
          *  time duration for steaming is 0.1 which is assigned for Molecular Dynamics and time duration for collision is ignored
          * 
          *  
          */


       for(i=0;i<=9999;i++)
 {
       double temp;
       double delta_t_MD=0.1; //time step

        temp=gsl_rng_uniform_pos(r);
        double rand_rx = (temp < 0.5) ? temp : ((-1) * temp ); // randomshift rx;

        temp=gsl_rng_uniform_pos(r);
        double rand_ry = (temp < 0.5) ? temp : ((-1) * temp ); // randomshift ry


       //steaming
      for(j=0;j<=3999;j++)
      {
            //printf("%d \n",j);
            printf("1 :%d \t  %f \t %f \n",j,solvent[j].rx,solvent[j].ry);
            solvent[j].rx=solvent[j].rx+solvent[j].vx * delta_t_MD;
            solvent[j].ry=solvent[j].ry+solvent[j].vy * delta_t_MD;

            printf("2: %d \t  %f \t %f \n",j,solvent[j].rx,solvent[j].ry);

    //randomshift
    solvent[j].rx=solvent[j].rx+rand_rx;
    solvent[j].ry=solvent[j].ry+rand_ry;
    printf("3:  %d \t  %f \t %f \n",j,solvent[j].rx,solvent[j].ry);
    // periodic boundary condition
    PBC(&solvent[j]);
    printf("4: %d \t  %f \t %f \n",j,solvent[j].rx,solvent[j].ry);
    solvent[j].cell=20*(ceil(solvent[j].ry)-1)+ceil(solvent[j].rx);
    printf("5: %d \t  %f \t %f \t %f \t %f \n",j,solvent[j].rx,solvent[j].ry,ceil(solvent[j].rx),ceil(solvent[j].ry));
    index = (int)(solvent[j].cell);
    //printf("%d \t %d \t %f \t %f \t %f \n",j,index,solvent[j].cell,ceil(solvent[j].rx),ceil(solvent[j].ry));

    if ((index>=0) &&(index<=400))
    {
    grid[index].vcm_x=grid[index].vcm_x+solvent[i].vx;
    grid[index].vcm_y=grid[index].vcm_y+solvent[i].vy;
    grid[index].number=grid[index].number+1;
   }
}


// caculating vcm

for (k=0;k<=399;k++)
{
    if (grid[k].number >= 1)
    {
        grid[k].vcm_x=grid[k].vcm_x/grid[k].number;
        grid[k].vcm_y=grid[k].vcm_y/grid[k].number; 
    }

    double temp;
    temp=gsl_rng_uniform_pos(r);
    grid[k].roation_angle = (temp < 0.5) ? alpha : ((-1) * alpha ); 
}



//collsion


 }

gsl_rng_free (r); // free RNG

return 0;

}

抱歉有点长，所以没有先放进去。有些东西还没有完成，但程序框架已经搭建好了。

我的程序是这样的:

    printf("4: %d \t  %f \t %f \n",j,solvent[j].rx,solvent[i].ry);
    solvent[j].cell=20*(floor(solvent[j].ry))+floor(solvent[j].rx)+1;
    printf("5: %d \t  %f \t %f \t %f \t %f \n",j,solvent[j].rx,solvent[i].ry,floor(solvent[j].rx),floor(solvent[j].ry));

尽管我想要更多的东西是错误的，下面是我选择输出的某些部分:

4: 3993       3.851240   17.047031 
5: 3993       3.851240   17.047031   3.000000    9.000000 

虽然floor(solvent[j].rx)是正确的，但floor(solvent[j].ry)是完全错误的。

我的程序最终结果是

 Segmentation fault (core dumped)


  ------------------
  (program exited with code: 139)

如何解决这个问题？我使用的东西有错吗？

为了进一步测试问题，我尝试了 ceil() 函数，程序和结果如下

程序:

    printf("4: %d \t  %f \t %f \n",j,solvent[j].rx,solvent[i].ry);
    solvent[j].cell=20*(ceil(solvent[j].ry)-1)+ceil(solvent[j].rx);
    printf("5: %d \t  %f \t %f \t %f \t %f \n",j,solvent[j].rx,solvent[i].ry,ceil(solvent[j].rx),ceil(solvent[j].ry));

结果:

2: 3999       14.571205      2.837654 
4: 3999       14.571205      2.837654 
5: 3999       14.571205      2.837654    15.000000   15.000000 

所以用最近的输出作为例子来说明我想要的结果是:

a=14.571205，ceil(a)=15.00 b= 2.837654 ， ceil(b)=3.00 不是输出中的 15.000

问题变得更糟的是，当我只使用 a 和 b 来测试 ceil() 时，一切看起来都很完美:

程序:

 #include <stdio.h>
 #include <math.h>

 int main()
{
        double a= 14.571205;  
        double b= 2.837654 ;  

        printf("%f \t %f \n",ceil(a),ceil(b));
        return 0;
 }

输出: 15.000000 3.000000

我在 Linux Ubuntu 中使用 GCC。

================================================== =================================

问题已解决。

真正致命的问题在这里

    if ((index>=0) &&(index<=400))
    {
    grid[index].vcm_x=grid[index].vcm_x+solvent[i].vx;
    grid[index].vcm_y=grid[index].vcm_y+solvent[i].vy;
    grid[index].number=grid[index].number+1;
   }
}

solvent[i].vy 和 solvent[i].vx 均应将 i 更改为 j。

正如 @Jon 和 @Blckknght @Eric Lippert 所指出的。

我显然陷入了错误的陷阱并误导了@ouah和@Blckknght，事实上，输出句子确实有问题，但它们不会导致prorame崩溃，只有超出边界才会这样做。

谢谢大家!

我确实想分享 Eric Lippert 的评论，它有力而富有洞察力:

More generally, the problem you have is called "select isn't broken" by experienced programmers. That is, you have done something completely wrong, you cannot figure out what you did wrong, and so you conclude that a piece of basic, easy infrastructure written by experts and tested extensively is wrong. I assure you, floor and ceil do exactly what they are supposed to do. They are not wrong. The problem lies somewhere in your code, not in the standard library. – Eric Lippert

Answer 1

您的数组声明为:

solvent_p solvent[4000];

但是你有这个循环:

 for(i=0;i<=9999;i++)

内部调用此函数:

 printf("1 :%d \t  %f \t %f \n",j,solvent[j].rx,solvent[i].ry);

这意味着您正在访问越界数组元素。

编辑:

OP 测试用例已被编辑以修复越界访问。

我的第二个建议是检查index值(用于索引grid数组)永远不会超过grid数组的元素数量分配后:index = (int)(solvent[j].cell)。