|
|
|
@ -2,7 +2,7 @@ |
|
|
|
#include <stdlib.h> |
|
|
|
#include <math.h> |
|
|
|
typedef struct { |
|
|
|
int N; |
|
|
|
double N; |
|
|
|
double d; |
|
|
|
double h; |
|
|
|
double xmin, xmax, tmin, tmax; |
|
|
|
@ -19,12 +19,13 @@ void ustart(soli_param *sp) { |
|
|
|
double jmax = sp->jmax; |
|
|
|
double h = sp->h; |
|
|
|
double x = sp->xmin; |
|
|
|
|
|
|
|
//FILE *fp2=fopen("cosh.csv","w"); |
|
|
|
for (int j = 0; j <= jmax; j++) { |
|
|
|
|
|
|
|
u0[j] = -N * (N + 1) / (pow(cosh(x), 2)); |
|
|
|
x = x + h; |
|
|
|
printf("%f\n", u0[j]); |
|
|
|
|
|
|
|
//fprintf(fp2,"%f\t%f\n", x,u0[j]); |
|
|
|
x = x + h; |
|
|
|
|
|
|
|
} |
|
|
|
|
|
|
|
@ -33,30 +34,32 @@ void ujeins(soli_param *sp) { |
|
|
|
double jmax = sp->jmax; |
|
|
|
double d = sp->d; |
|
|
|
double h = sp->h; |
|
|
|
for (int j = 0; j <= jmax; j++) { |
|
|
|
double x=sp->xmin; |
|
|
|
for (int j = 0; j < jmax; j++) { |
|
|
|
if (j == 1) { |
|
|
|
uj1[j] = (d * 6 * u0[j] * (u0[j + 1] - u0[j - 1]) / (2 * h)) |
|
|
|
- (d * (u0[j + 2] + 2 * u0[j - 1]) / (2 * pow(h, 3))) |
|
|
|
- (d * (u0[j + 2] + 2 * u0[j - 1]-2*u0[j+1]) / (2 * pow(h, 3))) |
|
|
|
+ u0[j]; |
|
|
|
|
|
|
|
} else if (j == 0) { |
|
|
|
uj1[j] = (d * 6 * u0[j] * (u0[j + 1]) / (2 * h)) |
|
|
|
- (d * (u0[j + 2]) / (2 * pow(h, 3))) + u0[j]; |
|
|
|
- (d * (u0[j + 2]-2*u0[j+1]) / (2 * pow(h, 3))) + u0[j]; |
|
|
|
|
|
|
|
} else if (j == jmax - 1) { |
|
|
|
} else if (j == jmax - 2) { |
|
|
|
uj1[j] = (d * 6 * u0[j] * (u0[j + 1] - u0[j - 1]) / (2 * h)) |
|
|
|
- (d * (2 * u0[j - 1] - u0[j - 2]) / (2 * pow(h, 3))) |
|
|
|
- (d * (-2*u0[j+1]+2 * u0[j - 1] - u0[j - 2]) / (2 * pow(h, 3))) |
|
|
|
+ u0[j]; |
|
|
|
|
|
|
|
} else if (j == jmax) { |
|
|
|
} else if (j == jmax-1) { |
|
|
|
uj1[j] = (d * 6 * u0[j] * (-u0[j - 1]) / (2 * h)) |
|
|
|
- (d * (2 * u0[j - 1] - u0[j - 2]) / (2 * pow(h, 3))) |
|
|
|
+ u0[j]; |
|
|
|
} else { |
|
|
|
uj1[j] = (d * 6 * u0[j] * (u0[j + 1] - u0[j - 1]) / (2 * h)) |
|
|
|
- (d * (u0[j + 2] + 2 * u0[j - 1] - u0[j - 2]) |
|
|
|
uj1[j] = ((d * 6 * u0[j] * (u0[j + 1] - u0[j - 1])) / (2 * h)) |
|
|
|
- (d * (u0[j + 2] -2*u0[j+1]+ 2 * u0[j - 1] - u0[j - 2]) |
|
|
|
/ (2 * pow(h, 3))) + u0[j]; |
|
|
|
} |
|
|
|
|
|
|
|
|
|
|
|
} |
|
|
|
|
|
|
|
@ -69,61 +72,76 @@ void nextu(soli_param *sp) { |
|
|
|
double * unm = u0; |
|
|
|
double * un = uj1; |
|
|
|
double unp[jmax]; |
|
|
|
for (int n = 0; n <= nmax; n++) { |
|
|
|
for (int j = 0; j <= jmax; j++) { |
|
|
|
double t=sp->tmin; |
|
|
|
double x=sp->xmin; |
|
|
|
int n; |
|
|
|
|
|
|
|
for ( n = 2; n < nmax; n++) { |
|
|
|
|
|
|
|
for (int j = 0; j < jmax; j++) { |
|
|
|
if (j == 1) { |
|
|
|
unp[j] = 2*d*( |
|
|
|
/*unp[j] = 2*d*( |
|
|
|
(un[j + 1] + un[j] + un[j - 1]) |
|
|
|
* (un[j + 1] - un[j - 1]) / h |
|
|
|
- (un[j + 2] - 2 * un[j + 1] |
|
|
|
+ 2 * un[j - 1]) |
|
|
|
/ (2 * pow(h, 3))) + unm[j]; |
|
|
|
/ (2 * pow(h, 3))) + unm[j];*/ |
|
|
|
unp[j]=(2*(d/h) *(un[j+1]+un[j]+un[j-1])*(un[j+1]-un[j-1]))-(d/((h*h*h))*(un[j+2]-2*un[j+1]+2*un[j-1]) ) +unm[j]; |
|
|
|
|
|
|
|
} else if (j == 0) { |
|
|
|
unp[j] = 2 |
|
|
|
/*unp[j] = 2 |
|
|
|
* d*( |
|
|
|
(un[j + 1] + un[j]) * (un[j + 1]) / h |
|
|
|
- (un[j + 2] - 2 * un[j + 1]) |
|
|
|
/ (2 * pow(h, 3))) + unm[j]; |
|
|
|
/ (2 * pow(h, 3))) + unm[j];*/ |
|
|
|
unp[j]=(2*(d/h) *(un[j+1]+un[j])*(un[j+1]))-(d/(h*h*h)*(un[j+2]-2*un[j+1]) ) +unm[j]; |
|
|
|
|
|
|
|
} else if (j == jmax - 1) { |
|
|
|
unp[j] = 2 |
|
|
|
} else if (j == jmax - 2) { |
|
|
|
/*unp[j] = 2 |
|
|
|
* d*( |
|
|
|
(un[j + 1] + un[j] + un[j - 1]) |
|
|
|
* (un[j + 1] - un[j - 1]) / h |
|
|
|
- (-2 * un[j + 1] + 2 * un[j - 1] |
|
|
|
- un[j - 2]) / (2 * pow(h, 3))) |
|
|
|
+ unm[j]; |
|
|
|
+ unm[j];*/ |
|
|
|
unp[j]=(2*(d/h) *(un[j+1]+un[j]+un[j-1])*(un[j+1]-un[j-1]))-(d/((h*h*h))*(-2*un[j+1]+2*un[j-1]-un[j-2]) ) +unm[j]; |
|
|
|
|
|
|
|
} else if (j == jmax) { |
|
|
|
unp[j] = 2 |
|
|
|
} else if (j == jmax - 1) { |
|
|
|
/*unp[j] = 2 |
|
|
|
* d*( |
|
|
|
(un[j] + un[j - 1]) * (-un[j - 1]) / h |
|
|
|
- (2 * un[j - 1] - un[j - 2]) |
|
|
|
/ (2 * pow(h, 3))) + unm[j]; |
|
|
|
/ (2 * pow(h, 3))) + unm[j];*/ |
|
|
|
unp[j]=(2*(d/h) *(un[j]+un[j-1])*(-un[j-1]))-(d/((h*h*h))*(2*un[j-1]-un[j-2]) ) +unm[j]; |
|
|
|
|
|
|
|
} else { |
|
|
|
unp[j] = 2 |
|
|
|
/*unp[j] = 2 |
|
|
|
* d*( |
|
|
|
(un[j + 1] + un[j] + un[j - 1]) |
|
|
|
* (un[j + 1] - un[j - 1]) / h |
|
|
|
- (un[j + 2] - 2 * un[j + 1] |
|
|
|
+ 2 * un[j - 1] - un[j - 2]) |
|
|
|
/ (2 * pow(h, 3))) + unm[j]; |
|
|
|
/ (2 * pow(h, 3))) + unm[j];*/ |
|
|
|
|
|
|
|
unp[j]=(2*(d/h) *(un[j+1]+un[j]+un[j-1])*(un[j+1]-un[j-1]))-(d/((h*h*h))*(un[j+2]-2*un[j+1]+2*un[j-1]-un[j-2]) ) +unm[j]; |
|
|
|
} |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
} |
|
|
|
|
|
|
|
unm=un; |
|
|
|
un=unp; |
|
|
|
x=sp->xmin; |
|
|
|
|
|
|
|
|
|
|
|
} |
|
|
|
FILE *fp= fopen("numplot2.csv", "w"); |
|
|
|
double x=sp->xmin; |
|
|
|
for(int j=0;j<=jmax;j++){ |
|
|
|
x=x+h; |
|
|
|
FILE *fp= fopen("numplot8.csv", "w"); |
|
|
|
printf("%f\t%f\n", d*n, d*nmax); |
|
|
|
for(int j=0;j<jmax;j++){ |
|
|
|
|
|
|
|
fprintf(fp, "%f\t%f\n",x,un[j]); |
|
|
|
printf( "%f\t%f\n",x,un[j]); |
|
|
|
//printf( "%f\t%f\n",x,un[j]); |
|
|
|
x=x+h; |
|
|
|
|
|
|
|
} |
|
|
|
fclose(fp); |
|
|
|
@ -131,20 +149,24 @@ void nextu(soli_param *sp) { |
|
|
|
|
|
|
|
|
|
|
|
int main(void) { |
|
|
|
soli_param sp; |
|
|
|
sp.N = 2; |
|
|
|
sp.h=0.1; |
|
|
|
sp.d=pow(sp.h,3)/3; |
|
|
|
sp.tmax=1; |
|
|
|
sp.tmin=0; |
|
|
|
sp.xmin=-30; |
|
|
|
sp.xmax=30; |
|
|
|
sp.jmax=(sp.xmax-sp.xmin)/sp.h; |
|
|
|
sp.nmax=(sp.tmax-sp.tmin)/sp.d; |
|
|
|
u0 = (double*) malloc(sp.jmax * sizeof(double)); |
|
|
|
ustart(&sp); |
|
|
|
uj1 = (double*) malloc(sp.jmax * sizeof(double)); |
|
|
|
ujeins(&sp); |
|
|
|
nextu(&sp); |
|
|
|
soli_param sp; |
|
|
|
sp.N = 1; |
|
|
|
sp.h=0.2; |
|
|
|
sp.d=pow(sp.h,3)/4; |
|
|
|
sp.tmax=1; |
|
|
|
sp.tmin=0; |
|
|
|
sp.xmin=-30; |
|
|
|
sp.xmax=30; |
|
|
|
sp.jmax=(sp.xmax-sp.xmin)/sp.h; |
|
|
|
sp.nmax=(sp.tmax-sp.tmin)/sp.d; |
|
|
|
//sp.nmax=500; |
|
|
|
printf("%d\n", sp.nmax); |
|
|
|
u0 = (double*) malloc(sp.jmax * sizeof(double)); |
|
|
|
ustart(&sp); |
|
|
|
uj1 = (double*) malloc(sp.jmax * sizeof(double)); |
|
|
|
ujeins(&sp); |
|
|
|
nextu(&sp); |
|
|
|
free(u0); |
|
|
|
free(uj1); |
|
|
|
|
|
|
|
} |