为什么文件输出的数据不对

huaijuliu

新手请教各位前辈。在init_f中试图输出两个向量f and hrhs，这两个向量都随着x的变动而变动，x从-4.5逐渐递增至1.5。两个文件ff1.dat and fhrhs1.dat的第一列均为x，第二列为对应的f and hrhs的值。由定义可知，f恒为零，而hrhs为抛物线曲线形状，但为什么fhrhs1.dat的第二列输出的hrhs的数据不对？请教各位。

1. 
   
2. #include <stdio.h>
   
3. #include <stdlib.h>
   
4. #include <math.h>
   
5. 
   
6. /* FMG Solver of the EHL circular contact */
   
7. #define pi 3.1415926535897931
   
8. typedef struct
   
9. {
   
10. double hx;
    
11. int ii;
    
12. double *p, *f;
    
13. double *hfi, *hrhs;
    
14. double *w;
    
15. double *K;
    
16. double *pjac;  /*old pressure for use in Jacobi relaxation*/
    
17. double rg;
    
18. double Hm;
    
19. } Level;
    
20. 
    
21. typedef struct
    
22. {
    
23. int nx0;
    
24. int maxlevel;
    
25. double xa,xb;
    
26. double h0;   /*global constant and work unit*/
    
27. Level *Lk;
    
28. } Stack;
    
29. 
    
30. void initialize(Stack *U, int nx0, int maxl, double xa, double xb, double h0)
    
31. {
    
32. /* initialize values in datastructure */
    
33. double hx;
    
34. Level *L;
    
35. int i,ii;
    
36. U->xa=xa;
    
37. U->xb=xb;
    
38. U->maxlevel=maxl;
    
39. U->h0=h0;
    
40. U->Lk=(Level *)calloc(maxl+1,sizeof(Level));
    
41. hx=(xb-xa)/nx0;
    
42. ii=nx0;
    
43. for (i=1; i<=maxl; i++)
    
44. {
    
45. L=U->Lk+i;
    
46. L->hx=hx;
    
47. L->ii=ii;
    
48. L->p =(double *)calloc(U->maxlevel+1,sizeof(double));
    
49. L->w =(double *)calloc(U->maxlevel+1,sizeof(double));
    
50. L->f =(double *)calloc(U->maxlevel+1,sizeof(double));
    
51. L->pjac =(double *)calloc(U->maxlevel+1,sizeof(double));
    
52. L->hfi =(double *)calloc(U->maxlevel+1,sizeof(double));
    
53. L->hrhs=(double *)calloc(U->maxlevel+1,sizeof(double));
    
54. L->K =(double *)calloc(U->maxlevel+1,sizeof(double));
    
55. printf("\n level: %2d ii=%4d, hx=%f",i,ii,hx);
    
56. hx*=0.5; ii*=2;
    
57. }
    
58. }
    
59. /********* SINGLE GRID ROUTINES *****/
    
60. void init_f(Stack *U, int lev)
    
61. {
    
62. int i;
    
63. Level *L;
    
64. double x;
    
65. L=U->Lk+lev;
    
66. for (i=0; i<=L->ii; i++)
    
67. {
    
68. x=U->xa+i*L->hx;
    
69. L->hrhs[i]=0.5*x*x;
    
70. L->f[i]=0.0;
    
71. }
    
72. printf("\n xa and hx are %f and %f",U->xa,L->hx);
    
73. FILE *ff1,*fhrhs1;
    
74. ff1=fopen("f1.dat","w");
    
75. fhrhs1=fopen("hrhs1.dat","w");
    
76. for (i=0; i<=L->ii; i++)
    
77. {
    
78. x=U->xa+i*L->hx;
    
79. fprintf(ff1,"%f %f\n",x,L->f[i]);
    
80. fprintf(fhrhs1,"%f %f\n",x,L->hrhs[i]);
    
81. }
    
82. fprintf (ff1, "\n");
    
83. fprintf (fhrhs1, "\n");
    
84. fclose(ff1);
    
85. fclose(fhrhs1);
    
86. L->rg=-pi/2.0;
    
87. }
    
88. 
    
89. /********** MAIN PROGRAM **********/
    
90. void main()
    
91. {
    
92. Stack U;
    
93. int i;
    
94. initialize(&U,128,1,-4.5,1.5,-0.5);//initialize(Stack *U, int nx0, int maxl, double xa, double xb, double h0)
    
95. init_f(&U,1);
    
96. }
    
97. 
    
98. 
    

复制代码

huaijuliu

问题已解决 感谢关注
initialize中的几条中的maxlevel改为ii
L->p =(double *)calloc(U->maxlevel+1,sizeof(double));

L->w =(double *)calloc(U->maxlevel+1,sizeof(double));

L->f =(double *)calloc(U->maxlevel+1,sizeof(double));

L->pjac =(double *)calloc(U->maxlevel+1,sizeof(double));

L->hfi =(double *)calloc(U->maxlevel+1,sizeof(double));

L->hrhs=(double *)calloc(U->maxlevel+1,sizeof(double));

L->K =(double *)calloc(U->maxlevel+1,sizeof(double));