Sphere Temperature in a Changing Bath: Can You Help Solve This PDE?

  • Thread starter Thread starter fluidistic
  • Start date Start date
  • Tags Tags
    Pde Sphere
Join the discussion
Registration is free. Ask a follow-up in this thread, or start your own.
21 replies · 4K views
fluidistic
Gold Member
Messages
3,934
Reaction score
286

Homework Statement


A sphere of radius R at temperature T=0 is put into a bath at time t=0 whose temperature is [itex]T_0[/itex].
Calculate the temperature inside the sphere [itex]\forall t \geq 0[/itex], [itex]T(\vec x ,t )[/itex].

Homework Equations


Heat equation: [tex]\frac{\partial T }{\partial t} \cdot \frac{1}{\kappa} -\triangle T =0[/tex]


The Attempt at a Solution


I will use separation of variables as well as making the assumption that due to the symmetry of the problem, T will only depend on r and t and not on theta and phi (I'm talking about spherical coordinates).
Thus [itex]T(\vec x , t ) = T(r,t)=\tau (t) R(r)[/itex].
Plugging this back into the PDE and taking the Laplacian in spherical coordinates, the PDE reduces to 2 ODE's, namely [tex]\begin{cases} \frac{1}{\kappa } \frac{\tau '}{\tau } =-C \\ \frac{2R'}{rR} + \frac{R''}{R}=-C \end{cases}[/tex] where C is a constant.
I solved the first ODE, the solution is [itex]\tau (t)=Ae^{-\kappa C t}[/itex].
I rewrote the second ODE into the form [itex]R''+\frac{2R'}{r}+CR=0[/itex]. To solve this DE I tried the Frobenius's method. Namely [itex]R(r)=\sum _{n=0}^ \infty a_n r^{n+\mu}[/itex]. The secular equation gave me [itex]\mu =-1[/itex] or [itex]\mu =0[/itex]. Since they differ by an integer I can only get one solution using this method and I'll have to use variation of parameters to get the linearly independent other solution.
So I have to take the lowest mu value, namely -1 here.
I reached that [itex]a_0[/itex] is arbitrary as well as [itex]a_1[/itex] but [itex]a_0[/itex] cannot be 0. Furthermore I obtained the following recurrence relation: [itex]a_n=\frac{-Ca_{n-2}}{n^2-n}[/itex], [itex]\forall n \geq 2[/itex].
Choosing [itex]a_0=1[/itex] and [itex]a_1=0[/itex], I sought to obtain the general form of [itex]a_n[/itex]. But I was not successful.
I reached that [itex]\forall n \geq 1[/itex], [tex]a_n =\begin{cases} 0 \text{ if n is odd} \\ \frac{ (-1)^{n/2}C^{n/2}}{(n^2-n)[(n-2)^2-(n-2)]...(2^2-2)} \text{ if n is even} \end{cases}[/tex].
I'm basically stuck at rewriting the denominator of [itex]a_n[/itex] when n is even. Can somebody help me?
This looks pretty awful!
 
Physics news on Phys.org
dextercioby said:
How did you implement the 2 boundary conditions, i.e. T(R,t) = T_0 and T(r,t_0) = 0 ?

I did not touch the boundary conditions yet, I'm looking for the general solution first and then I'll use the boundary condition+the initial condition.
I think they should be T(R,t)=T_0 and T(r<R,0)=0.Edit: Hmm I think I see what you mean. My solution for the [itex]\tau (t)[/itex] does not seem to make any sense, right?
So basically the method of separation of variables fails here?
 
dextercioby said:
Can you make a change of variable in your ODE, i./o. going to Frobenius ? Try u(r) = r R(r).

Ok thanks a lot. I will try this suggestion. By the way what does i./o. mean? :)
 
dextercioby said:
It should probably be i/o, a web shortage for "instead of", just like irrep stands for "irreducible representation".

Ah ok!
Well that's brillant, I get the ODE [itex]u''+Cu=0[/itex]. At first glance C will have to be negative, I'm going to proceed further in a few... I must leave for now.
 
If I call [itex]C=-k^2[/itex] the solution reduces to [itex]R(r)= \begin{cases} c_1 \frac{e^{kr}}{r}+ c_2 \frac{e^{-kr}}{r} \text{ for when } k \in \mathbb{R} \\ c_1 \frac{\cos (k_cr)}{r} + c_2 \frac{\sin (k_cr )}{r} \text { if k } \in \mathbb{C} \\ c_1+\frac{c_2}{r} \text{ if k =0} \end{cases}[/itex] where [itex]k_c[/itex] stands for the imaginary part of k.
Now of course I must dicard 2 of these solutions and find the constant(s) c1 and c2. Apparently k cannot be real because both c1 and c2 must be 0 to keep a finite solution at the origin.
If k is complex, then c1 must be worth 0. I think [itex]R(r)=c_2 \frac{\sin (k_cr )}{r}[/itex] could be possible.
I think I can safely discard the case k=0 because c2 would have to be worth 0 and by intuition the dependence of R(r) must not be a constant.
So I don't really know how to write the answer. I think it's [itex]R(r)=\frac{c_2 \sin ( k_c r)}{r}[/itex] but I don't know how to deal with [itex]k_c[/itex], it's not necessarily worth [itex]\sqrt C[/itex]. I'm not very confident of what I've done so far... what do you think?
 
haruspex said:
Not sure why you think C should be negative. Indeed, I'd be surprised if your k could have a real part, so I would expect an answer c sin(r√C)/r.

Sorry for my post 7, I was wrong. C indeed seems positive because an answer of the form c sin(r√C)/r assumes C is positive, for this solution is valid when k is a complex number and if it's purely imaginary then k^2 is a negative real number and thus C is positive.
So far I've reached the general form of the solution. It's [itex]T(r,t)=Ae^{-\kappa C t} c_1 \frac{\sin (r \sqrt C) }{r}=\frac{c_2e^{-\kappa Ct}\sin (r \sqrt C)}{r}[/itex].
So I have 2 unknowns (C and c2) and 1 boundary +1 initial condition.
But my problem is the following:
[itex]T(R,t)=T_0 \Rightarrow c_2=\frac{RT_0}{\sin (R \sqrt C )}[/itex]. Okay so far.
But the initial condition gives me troubles: [itex]T(r,0)=0 \Rightarrow \frac{c_2 \sin (r \sqrt C )}{r}=0 \Rightarrow c_2=0[/itex] (that I discard) or [itex]\sin (r \sqrt C )=0[/itex]. So that [itex]r\sqrt C =n \pi[/itex] with [itex]n \in \mathbb{N}[/itex]. This makes C not a constant anymore but a variable that depends on r...
Sigh. I really don't know what's wrong. :frown:

Edit: I've just checked out and [itex]R(r)=\frac{c_2 \sin (\sqrt C r ) }{r}[/itex] does satisfies the ODE R''+(2/r)R'+CR=0, which is a good sign.

Edit 2: When I take the limit when r tends to 0 (when I go to the center of the sphere), I get [itex]T(0,t)=\sqrt C c_2 e^{-\kappa Ct}[/itex] which basically shows that it cools down (instead of heating up!) with time and at t=0 is worth [itex]\sqrt C c_2[/itex] which must equal 0 by the initial condition. So something is very wrong, but I don't see what's wrong.
 
Last edited:
Clearly as t increases it will tend to be T0 everywhere. So there should be a constant term T0, and the general solution will be T(r,t) = T0-ƩAke-κk2tsin(kr)/r
Haven't tried to spot where you lost the T=constant solution to the ODE.
PS: I'm not suggesting k is an integer there. The set of k values is to be determined.
 
Note that if you add or subtract any constant to [itex]T(\vec x ,t )[/itex] it will still satisfy the heat equation. Thus, you might try letting [itex]\tau(\vec x ,t )[/itex] = [itex]T(\vec x ,t )[/itex] - ##T_0##. It will have the same general form of solution that you found for [itex]T(\vec x ,t )[/itex] exept the boundary and initial conditions will change. You should find that you can satisfy the boundary condition at r = R for [itex]\tau(\vec x ,t )[/itex] with an infinite number of different choices for your constant C. These values of C will be related to the k's in haruspex's expression for the solution.
 
Thanks for the help guys, though I'm still stuck.
haruspex said:
Clearly as t increases it will tend to be T0 everywhere. So there should be a constant term T0, and the general solution will be T(r,t) = T0-ƩAke-κk2tsin(kr)/r
Haven't tried to spot where you lost the T=constant solution to the ODE.
PS: I'm not suggesting k is an integer there. The set of k values is to be determined.
So [itex]T(r,t)=T_0 - \sum _k A_k e^{-\kappa k^2 t} \frac{\sin (kr )}{r}[/itex]. The initial condition gives me [itex]T(r,0)=0 \Rightarrow T_0 - \sum _k A_k \frac{\sin (kr )}{r}=0[/itex]. I don't know if it's possible that a constant can be written as a sum of a sinus function divided by the argument of the sinus function. Maybe this is a Fourier series, I'm not really sure. But I don't know how to get the values of k from it.
TSny said:
Note that if you add or subtract any constant to [itex]T(\vec x ,t )[/itex] it will still satisfy the heat equation. Thus, you might try letting [itex]\tau(\vec x ,t )[/itex] = [itex]T(\vec x ,t )[/itex] - ##T_0##. It will have the same general form of solution that you found for [itex]T(\vec x ,t )[/itex] exept the boundary and initial conditions will change. You should find that you can satisfy the boundary condition at r = R for [itex]\tau(\vec x ,t )[/itex] with an infinite number of different choices for your constant C. These values of C will be related to the k's in haruspex's expression for the solution.
Using [itex]\tilde T(\vec x ,t ) =T(\vec x , t ) + T_0[/itex], the boundary condition gives me [itex]\frac{c_2e^{-\kappa C t} \sin (R \sqrt C )}{R}=2T_0[/itex]. Which is impossible because the left side depends on t while the right side does not.
The initial condition gives me [itex]\frac{c_2 \sin ( r \sqrt C )}{r}-T_0=0[/itex] so that the constant is a function of r, impossible.

Yeah it's strange that my method did not reach the fact that adding a constant to the general solution would still make it a general solution.
 
fluidistic said:
So [itex]T(r,t)=T_0 - \sum _k A_k e^{-\kappa k^2 t} \frac{\sin (kr )}{r}[/itex]. The initial condition gives me [itex]T(r,0)=0 \Rightarrow T_0 - \sum _k A_k \frac{\sin (kr )}{r}=0[/itex]. I don't know if it's possible that a constant can be written as a sum of a sinus function divided by the argument of the sinus function. Maybe this is a Fourier series, I'm not really sure. But I don't know how to get the values of k from it.
We'll come back to that. Let's look at the other boundary condition first.
Using [itex]\tilde T(\vec x ,t ) =T(\vec x , t ) + T_0[/itex], the boundary condition gives me [itex]\frac{c_2e^{-\kappa C t} \sin (R \sqrt C )}{R}=2T_0[/itex].
No, r = R should give you:
[itex]T(R,t)=T_0 - \sum _k A_k e^{-\kappa k^2 t} \frac{\sin (kR )}{R} = T_0[/itex]
[itex]\sum _k A_k e^{-\kappa k^2 t} \frac{\sin (kR )}{R} = 0[/itex]
That should tell you what the k values are.
 
fluidistic said:
Using [itex]\tilde T(\vec x ,t ) =T(\vec x , t ) + T_0[/itex], the boundary condition gives me [itex]\frac{c_2e^{-\kappa C t} \sin (R \sqrt{C })}{R}=2T_0[/itex].

Try shifting the temperature the other way. Let

[itex]\tilde T(\vec x ,t ) =T(\vec x , t ) - T_0[/itex]

Then reconsider the boudary condition at r = R and find the possible values of [itex]\sqrt{C}[/itex] ## \equiv k## that will satisfy the boundary condition. [Whoops! I see haruspex already made essentially the same comment.]
 
Thanks once more seriously guys! I'm sorry for my careless arithmetics errors.
So I get [itex]k=\frac{n\pi}{R}[/itex] with n=0, 1, 2, ... This makes [itex]C=\left ( \frac{n\pi}{R} \right ) ^2[/itex].
The initial condition [itex]\tilde T (r,0)=0[/itex] gives me [itex]\sum _{n=0}^\infty A_n \frac{\sin \left ( \frac{n\pi}{R} \right ) ^2 }{r}=T_0[/itex]. So I must get those A_n coefficients and I'd be done.
The trick is to convert the infinite series into an integral? The limits of the integral would be 0 and R I guess and it will be with respect to r. I'm not 100% sure how to proceed there.
 
fluidistic said:
So I get [itex]k=\frac{n\pi}{R}[/itex] with n=0, 1, 2, ... This makes [itex]C=\left ( \frac{n\pi}{R} \right ) ^2[/itex].
Good. But you should be able to see that you can ignore the n = 0 term.
The initial condition [itex]\tilde T (r,0)=0[/itex] gives me [itex]\sum _{n=0}^\infty A_n \frac{\sin \left ( \frac{n\pi}{R} \right ) ^2 }{r}=T_0[/itex]. So I must get those A_n coefficients and I'd be done.
Remember, the argument of the sine function is ##\sqrt{C}r##.
The trick is to convert the infinite series into an integral? The limits of the integral would be 0 and R I guess and it will be with respect to r. I'm not 100% sure how to proceed there.
What type of infinite series do you get on the left if you multiply both sides of the equation [itex]\sum _{n=1}^\infty A_n ...=T_0[/itex] by ##r##?
 
TSny said:
Good. But you should be able to see that you can ignore the n = 0 term.
True, because we discarted the case k=0 already, right?
Remember, the argument of the sine function is ##\sqrt{C}r##.
Yes, I made a latex typo and forgot to include the "r" part.

What type of infinite series do you get on the left if you multiply both sides of the equation [itex]\sum _{n=1}^\infty A_n ...=T_0[/itex] by ##r##?
Hmm I don't really know. Looks like a sine Fourier series.
 
fluidistic said:
True, because we discarded the case k=0 already, right?

Right, k = 0 will give zero contribution due to sin(kr) = 0 when k = 0.
Yes, I made a latex typo and forgot to include the "r" part.
And also forgot to take the square root of C.
Hmm I don't really know. Looks like a sine Fourier series.

Yes. There's a standard way to find the coefficients of a Fourier series.
 
TSny said:
And also forgot to take the square root of C.
Oh right, another latex typo I did not see.
Yes. There's a standard way to find the coefficients of a Fourier series.
I look at http://mathworld.wolfram.com/FourierSineSeries.html. So [itex]A_n=\frac{2}{R} \int _0 ^R r T_0 \sin \left ( \frac{n\pi r}{R} \right ) dr[/itex].
I solved that integral by parts. I reach that [itex]A_n=(-1)^{n+1} \frac{2RT_0}{n\pi}[/itex].
Therefore the final answer of the problem would be... [tex]\tilde T(r,t)=T_0 \left [ 1+ \frac{2R}{\pi r} \sum _{n=1}^\infty \frac{(-1)^n}{n} \sin \left ( \frac{n\pi r}{R} \right ) \right ][/tex].
I hope I did not make any typo and more importantly, any mistake.
 
fluidistic said:
Therefore the final answer of the problem would be... [tex]\tilde T(r,t)=T_0 \left [ 1+ \frac{2R}{\pi r} \sum _{n=1}^\infty \frac{(-1)^n}{n} \sin \left ( \frac{n\pi r}{R} \right ) \right ][/tex].
I hope I did not make any typo and more importantly, any mistake.
Other than leaving out the e-kt term?
 
haruspex said:
Other than leaving out the e-kt term?

Oh right haha. :smile: My excuse is that it was almost 3 am.
Problem solved. Thank you very much guys.:approve: