[ff3d-users] how to use an external mesh

[Pierre Michaud]

Hi,
I have used ff3d for a few days to model heat conduction problems
(I am using the ff3d-2005-12-03 pre-compiled binary for Linux).
The examples attached are for a square plate of dimensions 
0.07x0.07x0.005 of uniform material.
I apply a uniform power density on the z=0 plane (Neumann condition) and 
a Robin BC on the z=0.005 side.
Lateral sides are adiabatic.
When I use a cartesian mesh, as in the example1.ff file attached, ff3d 
gives a correct result.
Actually, I have been able to model  more complex problems involving 
multiple layers
of different materials and non-uniform power density, and ff3d gives correct
results PROVIDED the mesh is fine enough.
However, ff3d uses a huge amount of memory. The finest cartesian mesh I 
could test is 160x160x30
and it uses up to 2.5 Gbytes (!) of virtual memory. I would like to test 
bigger meshes.
Using an external mesh may solve my problem because, for the problems I 
am interested in, I do not need
a uniform mesh. So I tried to use an external mesh generated with gmsh, 
but it seems that I am not using it correctly
(see example2.ff attached, it gives a wrong result). I also tried a mesh 
generated with tetgen, but I obtain the
same wrong result as with gmsh. Can someone tell me what I am doing wrong ?

Pierre


double width = 0.07;        // plate width
double thickness = 0.005;   // plate thickness
double k = 5;               // thermal conductivity
double h = 500;             // heat transfer coefficient
double P = 1;               // heat source power
double A = width^2;
double sourceside = width;  // heat source side length
double sourcearea = sourceside^2;
double R = (thickness/(k*A)) + 1/(h*A); // thermal resistance

vector n = (5,5,2);
vertex a = (-width/2,-width/2,0);
vertex b = (width/2,width/2,thickness);
mesh M = structured(n,a,b);

// for defining 2D power density (Robin BC)
function q = 
(P/sourcearea)*one((abs(x)<=(sourceside/2))&&(abs(y)<=(sourceside/2)));

solve(u) in M 
krylov(type=bicgstab),bicgstab(maxiter=1000000,epsilon=1e-5),memory(matrix=sparse)
{
  test(w)
  int(k*grad(u)*grad(w)) + int[M zmax](h*w*u) = int[M zmin](w*q);
};

cout << "   Exact temperature: " << R * P << "\n";
cout << "   Approximate temperature: " << u(0,0,0) << "\n";

example2.msh
Description: Mesh model

double width = 0.07;        // plate width
double thickness = 0.005;   // plate thickness
double k = 5;               // thermal conductivity
double h = 500;             // heat transfer coefficient
double P = 1;               // heat source power
double A = width^2;
double sourceside = width;  // heat source side length
double sourcearea = sourceside^2;
double R = (thickness/(k*A)) + 1/(h*A); // thermal resistance

mesh M = read(gmsh,"example2.msh");

// for defining 2D power density (Robin BC)
function q = 
(P/sourcearea)*one((abs(x)<=(sourceside/2))&&(abs(y)<=(sourceside/2)));

solve(u) in M 
krylov(type=bicgstab),bicgstab(maxiter=1000000,epsilon=1e-5),memory(matrix=sparse)
{
  test(w)
  // 16 and 14 are reference numbers for the z=0 and z=0.005 mesh boundaries
  int(k*grad(u)*grad(w)) + int[14](h*w*u) = int[16](w*q);
};

cout << "   Exact temperature: " << R * P << "\n";
cout << "   Approximate temperature: " << u(0,0,0) << "\n";