1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116
| \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\
#include <fstream>
#include <iostream>
#include <stdio.h>
#include <math.h>
#include "test.h"
#include <string>
#include <limits>
using namespace std;
double max (double a,double b)
{
if(a<b)
{return b;}
else
{return a;}
}
double min (double a,double b)
{
return -max(-a,-b);
}
double fplus(double a,double b,double c,double d)
{
double m1=max(a,-b);
m1=max(m1,0);
double m2=max(c,-d);
m2=max(m2,0);
return sqrt(m1*m1+ m2*m2);
}
double fmoins(double a,double b,double c,double d)
{
double m1=min(a,-b);
m1=min(m1,0);
double m2=min(c,-d);
m2=min(m2,0);
return sqrt(m1*m1+ m2*m2);
}
double delta (double x,double y, double z)
{return (x-y)/z;}
double u(double x, double y)
{ double res;
res = x*y;
return res;
}
double f(double x,double y, double z)
{return x*y + z;}
int main()
{
int Nx=10;
int Nt=10;
int Ny=10;
double X[10][10][10];
double stepx=1.0/Nx;
double stepy=1.0/Ny;
double stept=1.0/Nt;
double c;
FILE *projet;
projet=fopen("fin_projet.txt", "w");
for(int i=0;i<10;i++)
{
for(int j=0;j<10;j++)
{
X[0][i][j]= u(i*stepx,j*stepy);
}
}
for(int t=1;t<Nt;t++)
{
for (int i=1;i<Nx-1;i++)
{
for(int j=1;j<Ny-1;j++)
{
c= f(i*stepx, j*stepy, t*stept);
if (c >= 0)
{X[t][i][j]= X[t-1][i][j] + c* stept * fplus(delta(X[t-1][i+1][j],X[t-1][i][j],stepx),delta(X[t-1][i][j],X[t-1][i-1][j],stepx),delta(X[t-1][i][j+1],X[t-1][i][j],stepy),delta(X[t-1][i][j],X[t-1][i][j-1],stepy));}
else
{X[t][i][j]= X[t-1][i][j] + c* stept * fmoins(delta(X[t-1][i+1][j],X[t-1][i][j],stepx),delta(X[t-1][i][j],X[t-1][i-1][j],stepx),delta(X[t-1][i][j+1],X[t-1][i][j],stepy),delta(X[t-1][i][j],X[t-1][i][j-1],stepy));}
cout<<t<<i<<j<<" "<<X[t][i][j]<<endl;
fprintf(projet,"10.0%f 10.0%f 10.0%f \n",t,i,j,X[t][i][j]);
}
}
}
fclose(projet);
} |
Partager