2
$\begingroup$

I would like to perform a 3d FEM heat transfer calculation, which should have also included fluid dynamic simulation. But I want to simplify it into a 3d FEM without fluid dynamic simulation, but with a boundary condition which also requires to solve a 1d differential equation. A simplified example of "the simplified model" can be seen as follows.

Let's say I have an element and it looks like the following graphic enter image description here

region = 
 RegionDifference[Cuboid[{0, 0, 0}, {0.006536/2, 0.05, 0.0047}], 
  Cuboid[{0, 0, 0.0007}, {0.0015, 0.05, 0.0017}]];
DiscretizeRegion[region]

The equations is as follows

eq = With[{lambda = 1, c = 4200, rho = 1000, 
   w = 0.02}, {lambda Laplacian[u[x, y, z], {x, y, z}] == 
    NeumannValue[0, 
      x == 0 || x == 0.006536/2 || z == 0 || z == 0.05] + 
     NeumannValue[8 (u[x, y, z] - 20), y == 0] + 
     NeumannValue[0.8 (u[x, y, z] - 20), y == 0.0047] + 
     NeumannValue[
      1000 (u[x, y, z] - t[z]), (y == 0.0007 && 
         x < 0.0015) || (y == 1.7/1000 && 
         x < 0.0015) || (x == 0.0015 && (0.0007 < y < 0.0017))],
   -c rho 0.003 0.001 w D[t[z], z] == 1000 (t[z] - u[x, y, z]), 
   t[0] == 30}]

Basically, the outer surfaces have normal robin boundary conditions. The inner surface has also robin boundary conditions, but they require solving t[z]. This t[z] is only dependent on z.

Solving this equation with NDSolveValue

sol = NDSolveValue[eq, {u, t}, {x, y, z} \[Element] region, Method -> {"PDEDiscretization" -> {"FiniteElement", "MeshOptions" -> {"MaxCellMeasure" -> 10^(-6)}}}]

gives quite a lot of errors

Transpose::nmtx: The first two levels of {NDSolve`xs$20949,t} cannot be transposed.
Part::partw: Part 2 of Transpose[{NDSolve`xs$20949,t}] does not exist.
Transpose::nmtx: The first two levels of {NDSolve`xs$20949,Function[{x,y,z},30]} cannot be transposed.
Part::partw: Part 2 of Transpose[{NDSolve`xs$20949,Function[{x,y,z},30]}] does not exist.
Set::partw: Part 2 of Transpose[{NDSolve`xs$20949,t}] does not exist.
Rule::argr: Rule called with 1 argument; 2 arguments are expected.
Function::fpct: Too many parameters in {x,y,z} to be filled from Function[{x,y,z},30][z].
NDSolveValue::overdet: There are fewer dependent variables, {u[x,y,z]}, than equations, so the system is overdetermined.

Does anyone know how to solve this problem? A thousand thanks.

$\endgroup$
  • $\begingroup$ With the code NDSolve`ProcessEquations[... , DependentVariables-> {u,t}] , you have an unique and clear error message : "The function t[z] does not have the same number of arguments as independent variables (3)" ( NDSolve`ProcessEquations is the first step of the process of resolution of the pde used by NDSolve ) $\endgroup$ – andre314 Mar 22 '17 at 20:58
  • $\begingroup$ @andre so it means that NDSolve is not able to solve the couples 1d and 3D equations? $\endgroup$ – 407Peezy Mar 23 '17 at 0:07
  • $\begingroup$ Probably, it is not able (now) $\endgroup$ – andre314 Mar 23 '17 at 0:31
2
$\begingroup$

The idea is to transform t[z] to t[x,y,z]. Note that I also changed the predicate of the y==0.0047 to z==0.0047.

region = RegionDifference[
   Cuboid[{0, 0, 0}, {0.006536/2, 0.05, 0.0047}], 
   Cuboid[{0, 0, 0.0007}, {0.0015, 0.05, 0.0017}]];
eq = With[{lambda = 1, c = 4200, rho = 1000, 
    w = 0.02}, {lambda Laplacian[u[x, y, z], {x, y, z}] == 
     NeumannValue[0, 
       x == 0 || x == 0.006536/2 || z == 0 || z == 0.05] + 
      NeumannValue[8 (u[x, y, z] - 20), y == 0] + 
      NeumannValue[0.8 (u[x, y, z] - 20), z == 0.0047] + 
      NeumannValue[
       1000 (u[x, y, z] - t[x, y, z]), (y == 0.0007 && 
          x < 0.0015) || (y == 1.7/1000 && 
          x < 0.0015) || (x == 
           0.0015 && (0.0007 < y < 0.0017))], -c rho 0.003 0.001 w D[
       t[x, y, z], z] == 1000 (t[x, y, z] - u[x, y, z]), 
    t[x, y, 0] == 30}];
sol = NDSolveValue[eq, {u, t}, {x, y, z} \[Element] region, 
  Method -> {"PDEDiscretization" -> {"FiniteElement", 
      "MeshOptions" -> {"MaxCellMeasure" -> 10^(-6)}}}]

This gives a message:

NDSolveValue::femibcnd: No DirichletCondition or Robin-type NeumannValue was specified for {u}; the result is not unique up to a constant.

but returns two interpolating functions. See if this helps you.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.