Can Matlab give me an analytic solution for transient Heat diffusion?
조회 수: 5 (최근 30일)
이전 댓글 표시
I want an analytic solution to the transient Heat Diffusion equation in cylindrical coordinates. My equation is
dT/dt = d^2T/dr^2 +1/r dT/dr + d^2T/dz^2
This problem is cylindrically symmetric, so the theta term is zero.
Is there a way to solve this in MatLab? My prof tells me MatLab should be able to do this with the symbolic toolkit, but all I can find on it concerns ODEs.
The boundary conditions for this problem are T(0,z,t) = g(z), T(r,z,0)=0
댓글 수: 2
Shaihroz Khan
2015년 11월 30일
Listen I don't know but u should better try to visualize your problem in fluent rather than wasting your time matlab, coz for that u have to program the code and you wont find it here.
Im also having the same problem but in momentum transfer for getting velocity profile in a cylinder with velocity being a function of r and z and time.
Walter Roberson
2015년 11월 30일
Shaihroz, this was from two years ago, and discussion stalled because the poster was unable to define the desired behaviour along the axis of the cylinder. If a behaviour had been defined we might have been able to program it.
답변 (1개)
Walter Roberson
2013년 11월 21일
Older versions of MATLAB used the Maple symbolic package, which has a PDE solver. I have not encountered a PDE solver for the MuPAD symbolic engine as yet, as partly because I do not have that Toolbox.
Maple believes that the only solution is if g(z) = 0, with T(r,z,t) = 0.
Maple is able to do a partial resolution if you omit the T(0,z,t)=g(z) initial condition. That is partly because solving PDE typically depends critically on details of all of the functions, but it is also because it encounters a singularity at r=0. The formula has a 1/r which at r=0 would give "undefined" leaving the rest of the formula undefined. When that first initial condition is omitted,
T(r, z, t) = -_C5 * _C1 * (_C3 * exp(_c[2]^(1/2) * z - t * RootOf(BesselJ(0, _Z * r) * BesselY(1, _Z * r) - BesselJ(1, _Z * r) * BesselY(0, _Z * r))^2 + t * _c[2]) + _C4 * exp(-_c[2]^(1/2) * z - t * RootOf(BesselJ(0, _Z * r) * BesselY(1, _Z * r) - BesselJ(1, _Z * r) * BesselY(0, _Z * r))^2 + t * _c[2])) * (-BesselY(1, RootOf(BesselJ(0, _Z * r) * BesselY(1, _Z * r) - BesselJ(1, _Z * r) * BesselY(0, _Z * r)) * r) * BesselJ(0, csgn(RootOf(BesselJ(0, _Z * r) * BesselY(1, _Z*r) - BesselJ(1, _Z*r) * BesselY(0, _Z*r))) * RootOf(BesselJ(0, _Z*r) * BesselY(1, _Z * r) - BesselJ(1, _Z*r) * BesselY(0, _Z*r)) * r) + BesselJ(1, RootOf(BesselJ(0, _Z*r) * BesselY(1, _Z * r) - BesselJ(1, _Z*r) * BesselY(0, _Z*r)) * r) * BesselY(0, csgn(RootOf(BesselJ(0, _Z*r) * BesselY(1, _Z * r) - BesselJ(1, _Z*r) * BesselY(0, _Z*r))) * RootOf(BesselJ(0, _Z*r) * BesselY(1, _Z*r) - BesselJ(1, _Z * r) * BesselY(0, _Z*r)) * r)) / BesselY(1, RootOf(BesselJ(0, _Z*r) * BesselY(1, _Z*r) - BesselJ(1, _Z*r) * BesselY(0, _Z*r)) * r)
This is Maple notation.
_C1, _c[2], _C3, _C4, and _C5 should all be read as the names of constants of integration that need to be resolved according to further initial conditions.
MuPAD uses lower-case names for the Bessel functions, such as bessely where Maple uses BesselY
The form RootOf(f(_Z)) in Maple translates into MuPAD as RootOf(f(z1), z1) . It stands for the set of values of the variable _Z such that f(_Z) is 0 -- the "roots" of the expression.
At r=0, BesselY(1, _Z*r) would be BesselY(1, 0) which is undefined. lim x => -0 BesselY(1, x) is +infinity and lim x <= +0 BesselY(1,x) is -infinity so there is an unresolvable singularity.
Perhaps T(0,z,t) = g(z) is not really a "boundary condition" but rather a definition overriding the behavior at r = 0 ?
댓글 수: 2
Walter Roberson
2013년 11월 22일
Suppose you define a function
y = 6/x
Then as you approach x = 0 from above, y goes to +infinity, but as you approach x = 0 from below, y goes to -infinity. So in this example function, y is well defined for every real number except x = 0 exactly, where it must be undefined, as it has a fundamental singularity there.
As opposed to the example
y = 6 / abs(x)
then as x approaches 0 from above , y goes to +infinity, an as x approaches 0 from below, y goes the +infinity. Right at x = 0 it has a singularity, but it is a "removable singularity", because if you were to _define y(0) as +infinity then there would be continuity as you flowed upward toward y = infinity as you approach x = 0 from below, then you stay at infinity at x = 0 exactly, then you come down from infinity as you get further positive away from x = 0.
Thus, if what you were to do would be to define your function as:
t(r,z,t) = piecewise( r == 0, -infinity, otherwise dT/dt = d^2T/dr^2 +1/r dT/dr + d^2T/dz^2 )
then the function would be well-defined at every point because at what would otherwise be the plane of singularity at r = 0 you have defined a result.
The difficulty with doing this is that you would still have a singularity at r = 0 because although the function is well defined for every r = -delta_R, delta_R > 0, and is well defined for every r = 0, and is well defined for every r = +delta_R, delta_R = -infinity, there is no continuity as you go from negative r to r = 0: you go from +infinity to -infinity.
What you can do is to say that the domain of the function is r >= 0, which causes the problems with negative r to vanish because the function would just not apply, just as it is not meaningful to "take the square root of a walnut". Then you are still stuck with what to do at r = 0 exactly. If you define T(0,z,t) = g(z) then you define a value for each (r,z,t), r >= 0, but unfortunately unless g(z) = -infinity always, then you still get a discontinuity because T(r,z,t) goes to -infinity as r approaches 0 from above and hence any finite value defined at r = 0 exactly would not be adjacent to the -infinity.
You need to rethink what should happen at r = 0 exactly. If r indicates a distance from the central axis of the cylinder, then g(z) is dependent only on the height and not on the time: is that truly reasonable? The value for every molecule touching those in the central core has a value that is dependent on height and time both; if any kind of solid is being modeled, one would expect those time-dependent temperatures to affect temperatures at the pole. Or is there a theorem that states that the inflow of heat from the adjacent molecules will always exactly balance the outflow of heat?
참고 항목
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!