菜鸟请各位高手帮忙了!!谢谢
#0 0x08048cd5 in SPH_accelCalculation (numArray=Cannot access memory at address 0xbfe18af0) at SPH_accel.c:99
Cannot access memory at address 0xbfe18aec
这是通过gdb显示的段错误,但是不知道怎么会发生这种错误。程序见下:
void SPH_accel( int numArray, SPH_t *sph )
{
SPH_accelCalculation( sph->map->numArray, sph->node, sph->mater, sph->value, sph->pressure, sph->cell );
}
/*===========================================================================*/
void SPH_accelCalculation( int numArray, Node_t *node, Mater_t *mater, Value_t *value, Pressure_t *pressure, Cell_t *cell)
{
Value_t *IVal =value+cell->scopeList[numArray];
double *Iden= cell->density_c + numArray;
double *Ivel= cell->velocity_c + numArray;
double *IvelS= cell->velocityS_c + numArray;
double *IaWaX= cell->aWaX_c + numArray;
Mater_t *IMat = mater;
Pressure_t *Ipress = pressure;
double *Jden;
double *Jvel;
Mater_t *JMat;
double *Ju;
int *Jtype;
Pressure_t *Jpress;
int a, i;
double sxx, sxy, sxz, syy, syz, szz, Itxx, Itxy, Itxz, Ityy, Ityz, Itzz, Jtxx, Jtxy, Jtxz, Jtyy, Jtyz, Jtzz;
double hvcc, Ivcc, Jvcc, tdsdt, rhoij, p0[DIM], p[DIM], pp[DIM];
/*******: Initialize accelation *************/
for( i=0; i<DIM; ++i) IVal->accel[i]=0.0;
if(node->type==0){
for( a=0; a< cell->sumScope; a++){
Ju = cell->u_c+a;
if( Ju[0] > 0.0) {
Jden = cell->density_c+a;
Jvel = cell->velocity_c+a;
Jtype = cell->type_c+a;
JMat = mater+Jtype[0];
Jpress = pressure;
sxx=sxy=sxz=syy=syz=szz=Itxx=Itxy=Itxz=Ityy=Ityz=Itzz=Jtxx=Jtxy=Jtxz=Jtyy=Jtyz=Jtzz=0.0;
hvcc=Ivcc=Jvcc=tdsdt=rhoij=0.0;
p0[i]= 0.0;
p[i]=0.0;
pp[i]=0.0;
/*calcutlate SPH sum for shear tensor Tab = va,b +vb,a - 2/3delta_ab vc,c */
sxx = 2.0*(Ivel[0]-Jvel[0])*IaWaX[0]-(2.0/3.0*(Ivel[1]-Jvel[1])*IaWaX[1])-(2.0/3.0*(Ivel[2]-Jvel[2])*IaWaX[2]);
sxy = ( Ivel[0] - Jvel[0] ) * IaWaX[1] + ( Ivel[1] - Jvel[1] ) * IaWaX[0];
sxz = ( Ivel[0] - Jvel[0] ) * IaWaX[2] + ( Ivel[2] - Jvel[2] ) * IaWaX[0];
syy = 2.0*(Ivel[1]-Jvel[1])*IaWaX[1]-(2.0/3.0*(Ivel[0]-Jvel[0])*IaWaX[0])-(2.0/3.0*(Ivel[2]-Jvel[2])*IaWaX[2]);
syz = ( Ivel[1] - Jvel[1] ) * IaWaX[2] + ( Ivel[2] - Jvel[2] ) * IaWaX[1];
szz = 2.0*(Ivel[2]-Jvel[2])*IaWaX[2]-(2.0/3.0*(Ivel[1]-Jvel[1])*IaWaX[1])-(2.0/3.0*(Ivel[0]-Jvel[0])*IaWaX[0]);
Itxx += JMat->mass * sxx / Jden[i];
Jtxx += IMat->mass * sxx / Iden[i];
Itxy += JMat->mass * sxy / Jden[i];
Jtxy += IMat->mass * sxy / Iden[i];
Itxz += JMat->mass * sxz / Jden[i];
Jtxz += IMat->mass * sxz / Iden[i];
Ityy += JMat->mass * syy / Jden[i];
Jtyy += IMat->mass * syy / Iden[i];
Ityz += JMat->mass * syz / Jden[i];
Jtyz += IMat->mass * syz / Iden[i];
Itzz += JMat->mass * szz / Jden[i];
Jtzz += IMat->mass * szz / Iden[i];
/* calculate SPH sum for vc,c=dvx/dx + dvy/dy + dvz/dz */
hvcc += (Ivel[i]-Jvel[i]) * IaWaX[i];
Ivcc += JMat->mass * hvcc/Jden[i];
Jvcc += IMat->mass * hvcc/Iden[i];
/* calculate SPH sum for pressure force -p,a/density and viscous force (visco Tab), b/density */
rhoij = 1.0/(Iden[i]+Jden[i]);
//concered with J paticle
p0[i] = -(Ipress->p + Jpress->p)*IaWaX[i];
p[0] = p0[0]+((mater->visco*Itxy+mater->visco*Jtxy)*IaWaX[1]+(mater->visco*Itxz+mater->visco*Jtxz)*IaWaX[2]);
p[1] = p0[1]+((mater->visco*Itxy+mater->visco*Jtxy)*IaWaX[0]+(mater->visco*Ityz+mater->visco*Jtyz)*IaWaX[2]);
p[2] = p0[2]+((mater->visco*Itxz+mater->visco*Jtxz)*IaWaX[0]+(mater->visco*Ityz+mater->visco*Jtyz)*IaWaX[1]);
pp[i] = p[i] * rhoij;
printf("00000");
//Calculate accel for each direction
IVal->accel[i]+=(-JMat->mass)*pp[i]-IMat->beta*(Ivel[i]-IvelS[i])/(mater->density*IMat->epsilon);
//Consideration of gravity
IVal->accel[2]= IVal->accel[2]-JMat->gravity;
printf("a=%d, Iaccel[0]=%.19e, Iaccel[1]=%.19e, Iaccel[2]=%.19e\n", a, IVal->accel[0], IVal->accel[1], IVal->accel[2]);
printf("IaWaX[2]=%.19e, IaWaX[1]=%.19e, IaWaX[0]=%.19e \n", IaWaX[2], IaWaX[1], IaWaX[0]);
}
}
}
}