20 |
// A real-world test case from Stephane Del Pino's PhD Thesis, |
// A real-world test case from Stephane Del Pino's PhD Thesis, |
21 |
// paragraph 4.2.4: natural convection. |
// paragraph 4.2.4: natural convection. |
22 |
|
|
23 |
vector n = (20,20,20); |
vector n = (10,10,10); |
24 |
vertex A = (0,0,0); |
vertex A = (0,0,0); |
25 |
vertex B = (1,1,1); |
vertex B = (1,1,1); |
26 |
|
|
41 |
double K = 0.1; |
double K = 0.1; |
42 |
double eps = 1E-3; |
double eps = 1E-3; |
43 |
|
|
44 |
function Ki = 1-one(<1,0,0>); |
double i; |
45 |
|
for (i = 1;i <= 10;i = i+1) { |
46 |
|
|
47 |
|
cout << "========== Natural convection step " << i << "\n"; |
48 |
|
|
|
double i=0; |
|
|
for (i=0; i<10; i=i+1) { |
|
49 |
solve(theta) in Omega by M method(type=penalty), |
solve(theta) in Omega by M method(type=penalty), |
50 |
krylov(precond=diagonal), memory(matrix=sparse), cg(epsilon=1E-10) |
krylov(precond=diagonal), memory(matrix=sparse), cg(epsilon=1E-10) |
51 |
{ |
{ |
66 |
cg(epsilon=1E-5), memory(matrix=none) |
cg(epsilon=1E-5), memory(matrix=none) |
67 |
{ |
{ |
68 |
pde(u1) |
pde(u1) |
69 |
(1-Ki)*u1-div(grad(u1))+dx(p)=0; |
-div(grad(u1))+dx(p)=0; |
70 |
|
|
71 |
u1 = 0 on M; |
u1 = 0 on M; |
72 |
u1 = 0 on T; |
u1 = 0 on T; |
73 |
|
|
74 |
pde(u2) |
pde(u2) |
75 |
(1-Ki)*u2-div(grad(u2))+dy(p)=K*theta0; |
u2-div(grad(u2))+dy(p)=K*theta0; |
76 |
|
|
77 |
u2 = 0 on M; |
u2 = 0 on M; |
78 |
u2 = 0 on T; |
u2 = 0 on T; |
79 |
|
|
80 |
pde(u3) |
pde(u3) |
81 |
(1-Ki)*u3-div(grad(u3))+dz(p)=0; |
u3-div(grad(u3))+dz(p)=0; |
82 |
|
|
83 |
u3 = 0 on M; |
u3 = 0 on M; |
84 |
u3 = 0 on T; |
u3 = 0 on T; |
86 |
pde(p) |
pde(p) |
87 |
-div(eps*grad(p))+eps*p+dx(u1)+dy(u2)+dz(u3)=0; |
-div(eps*grad(p))+eps*p+dx(u1)+dy(u2)+dz(u3)=0; |
88 |
}; |
}; |
|
|
|
|
cout << "step " << i << " done\n"; |
|
|
save(medit,"natural-convection-u-".i, u1, M); |
|
|
save(medit,"natural-convection-v-".i, u2, M); |
|
|
save(medit,"natural-convection-w-".i, u3, M); |
|
|
save(medit,"natural-convection-p-".i, p, M); |
|
|
save(medit,"natural-convection-theta-".i, theta0, M); |
|
89 |
} |
} |
90 |
|
|
91 |
save(medit,"natural-convection-u-".i,M); |
save(medit,"natural-convection-u",M); |
92 |
save(medit,"natural-convection-v-".i,M); |
save(medit,"natural-convection-u", u1, M); |
93 |
save(medit,"natural-convection-w-".i,M); |
save(medit,"natural-convection-v",M); |
94 |
save(medit,"natural-convection-p-".i,M); |
save(medit,"natural-convection-v", u2, M); |
95 |
save(medit,"natural-convection-theta-".i,M); |
save(medit,"natural-convection-w",M); |
96 |
|
save(medit,"natural-convection-w", u3, M); |
97 |
|
save(medit,"natural-convection-p",M); |
98 |
|
save(medit,"natural-convection-p", p, M); |
99 |
|
save(medit,"natural-convection-theta",M); |
100 |
|
save(medit,"natural-convection-theta", theta0, M); |