The 1D thermal diffusion problem is as follows:
div[k grad(theta)] = ro * c * du
where: ro = density, k = coefficient of thermal conduction , c = specific heat at constant deformation, du = velocity.
Can anyone please tell me the example boundary condition for this problem, so that I can implement it in my FEM code?
Dear Abdullah,
in your equation, you most likely forgot the temperature gradient on the right hand side, i.e.
div[k grad(theta)] = rou * c * du * grat(theta)
If this is the case, then when you introduce the Peclet number, Pe=du * L / k * rho * c,
you can write, for example in 1D for theta=theta(x):
Pe*d(theta)/dx=d^2(theta)/dx^2
if you solve for x=(0,1) and use boundary conditions theta(0)=1 and theta(1)=0, then the analytical solution is
theta(x)=1-(1-Exp(Pe*x))/(1-Exp(Pe))
If you set k=x, then
Pe*d(theta)/dx=d/dx(x*d(theta)/dx)
In this case, for the same boundary conditions, the analytical solution is
theta(x)=1-x^Pe
Hope this helps,
Jure.
Hi, Abdullah !
Boundary conditions must correspond to the physical problem to be solved . If the above equation describes the process of stationary heat transfer , it is desirable to clarify .
1. The distribution function of power density of heat sources (the right side of the equation ) depends on the speed , but the speed of what?
2 . In the left and right boundary points in the simplest case it can be assigned one of three possible conditions - thermal insulation (adiabatic condition), certain fixed temperature value or condition of convective heat transfer. Selecting the most likely of these conditions is based on the physical characteristics of the process (adopted assumptions).
One option is make everything zero at boundaries and initials and then include a heat source in the equation ( e.g. gaussian source). or
at x=0, theta = 0 ( tempertature at surface of the earth )
d(theta)/dx at x=0 , is 30 degree/km
etc.
That depends on which BC you want to apply, i.e Isothermal BC or adiabatic BC.
1) I agree with Jure, you probably want ro*c*du*grad(theta) on the right hand side. Then, the left-hand side describes diffusion and the left-hand side describes convection (http://goo.gl/aGH6Ln)
2) If "k" and "du" are functions of "x", then you need the following form:
div [a*grad(theta)] = div [du * theta]
where "a" is the thermal diffusivity, a = k/(ro*c)
In this case, you can't just take the "k" out of the div operator. However, when you are testing your FEM implementation, use constant "a" and "du" to compare your FEM solution to an analytic solution. If "a" and "du" are variable, it may be impossible to find an analytic solution - that's why you need FEM!
3) If you're testing your FEM, I would suggest using a fixed temperature BC on one side (e.g. at x = 0) and an insulating BC on the other, because the way you implement natural and essential BCs is different (http://goo.gl/I3jlDx). That way, you can see whether you're implementing the two types of BCs correctly.
Here is an example problem:
At steady state, in one dimension only, with constant "a" and "du" :
d^2 (theta) / dx^2 - du/a * d(theta)/dx = 0
BCs: theta = theta0 at x=0 (essential BC)
d(theta)/dx = q at x = L (natural BC)
Solution: theta = theta0 + (q*a/du) * exp(-du*L/a) * ( exp(du*x/a) - 1 )
Compare the analytic solution to the FEM results to see if you've implemented it correctly.
Tobias Louw .. Thank you .. I have implemented time dependent Dirichlet boundary conditions on both sides of the domain as u(0)=u(L)=0 and Initial temperature at the middle of domain and its working correctly.
Thank You all for responding to the post, your answers always leads to increase in knowledge and helpful.
There's some useful examples together with Fortran source code in the 5th edition of "Programming the Finite Element Method". See http://tinyurl.com/o27v6s2 for a link to a description of the book. Perhaps your University library can source a copy for you?
- Badges
- Science topic