Lagrangian Mechanics and The Heat Equation

I have been noticing a lot of discussion on Lucia’s blog about an energy Balance Climate Model.

Two Box Model: Algebra for ver. 1 test.
Two Box Model: Now assuming surface temperature are ‘mixed’ values.
Two Box Model: Rough idea how to obtain parameters.
Two Box Models & The 2nd Law of Thermodynamics

Which started based on a post by Tamino, called “Two Box”

The biggest criticism I’ve scene about these models is that there is very little physics which is used to form these models and consequently they won’t capture much of the nonlinear dynamics.

The basic premises of the model is that the time lags are associated with the heat capacity of the ocean. Consequently, the ocean and some of the atmosphere is partitioned into boxes and the model tries to determine how these boxes are coupled.

The end result though is the time legs end up being associated with eignvectors which, project onto both boxes. Consequently the modes of the system do not coincide with the boxes. Moreover, it is desirable to minimize the number of model parameters which need to be fit and therefore principles of physics, should be used to derive a more realistic “Two box model” or should I say more generally “Two mode model” as the modes do not coincide with the boxes.

The most simple equation with regards to the transfer of energy is the heat equation:

1) $\frac{\partial u}{\partial t} -k\left(\frac{\partial^2u}{\partial x^2}+\frac{\partial^2u}{\partial y^2}+\frac{\partial^2u}{\partial z^2}\right)=0$

The heat equation is essentially a diffusion equation based on Brownian Motion. But we can adjust the constant $K$ to try and account for other heat exchange processes.

In one dimension the heat equation can be written as:

2) $\frac{\partial u}{\partial t} -k\frac{\partial^2u}{\partial y^2}=0$

In this equation the time derivative only depends on the spatial derivatives. Keeping a mind on the dynamics, the higher the frequency that the climate forcing is, the more shallow it will penetrate into the ocean because the more ocean you include the greater the thermal inertia (heat capacity) and consequently the greater the damping. Thus a good fit for the temperature response at a given frequency as a function of depth, may be an exponential function of the depth.

Thus, from this observation and because the mathematics are fairly simple consider a set of basis functions which the special component is exponentially decreasing from the surface. That is let:

3) $u=\sum a_i exp(-\lambda_i y)$

Plugging this into the heat equation one gets:

4) $\frac{du}{dt}=\sum \frac{d ( \ a_i)}{dt} exp(-\lambda_i )=k \lambda_i^2\sum a_i exp(-\lambda_i y)$

Where:

The constant $k(y)$ represents how quickly heat diffuses in the ocean at a given depth. Equating terms on the heat equation gives:

5) $k\lambda_i^2a_i-\frac{d (a_i)}{dt}=0$

Now let there be a function $c(y)$ , and let the quanity being diffused in the heat equation be the temperature. Then the energy constraint is given by:

6) $E=\int_0^{d} c(y)\sum a_i exp(-\lambda_i y)dy$

This constraint makes one of the generalized coordinates $a_i$ redundant.

Differentiating (5) with respect to time, gives:

7) $\frac{dE}{dt}=\int_0^{d} \sum c(y)\frac{d(a_i)}{dt} exp(-\lambda_i y)dy$

$\frac{dE}{dt}$ is determined by the forcing. Now some notes, with regards to the coordinate system:

1)I suggest taking the as the redundant basis function, the one that has the smallest value $\lambda$ .
2)Since this model is a linearization (or approximation), I suggest using as the temperature variable at a given depth, the amount by which it exceeds the average temperature at that depth. That is model the temperature anomaly instead of the actual temperature.
3)Given two then we are looking at the change in the heat transfer from the mean instead of the actual heat transfer.
4)These are only suggestions.

Now the Lagrangian, is essentially a cost function. In our cost function, we want to minimize the error in the derivatives (the square of equation (5) integrated over y) and as an additional constraint, I want to minimize the change in entropy with respect to y.

When a parcel of water rises because it is hotter then the surrounding and therefore less dense it will expand with constant entropy (adiabatic expansion) if it exchanges no heat with the surroundings. Therefore, the heat induced ocean currents will act in order to try and minimize the entropy change between adjacent regions. (Note the actual effect of this process is to increase entropy).

From thermodynamics:

8 ) $dU=T\,dS-p\,dV$

Therefore:

9) $dS=(dU-p\,dV)/T=\frac{1}{T}(\frac{1}{C(y)}dT-p(y) \frac{dV}{dP}dP)$

from which we can get:

10) $\frac{dS}{dy}=\frac{1}{T}(\frac{1}{C(y)}\frac{dT}{dy}-p(y) \frac{dV}{dP}\frac{dP}{dy})$

Substituting into (10) the expression for temperature (equation (3)) equation (10) gives:

11) $\frac{dS}{dy}=-\frac{1}{T}(\frac{1}{C(y)}\sum a_i\frac{-1}{\lambda_i} exp(\lambda_i y)+p(y) \frac{dV}{dP}\frac{dP}{dy})$

The entropy cost function is:

12) $\int_o^d (\frac{dS}{dy})^2dy=\int_0^d\frac{1}{T^2}(\frac{1}{C(y)}\sum a_i\frac{1}{\lambda_i} exp(\lambda_i y)+p(y) \frac{dV}{dP}\frac{dP}{dy})^2dy$

The Lagrangian is given by:

13) $\mathcal L(a_1,...a_n, ...;\dot a_1,..,\dot a_n, ...;t)\,.= \alpha \int_o^d (\frac{dS}{dy})^2dy+\sum(k\lambda_i^2a_i-\frac{d (a_i)}{dt})^2$

Where $\alpha$ is a stiffness parameter for the entropy.

Equation 13) is subject to the constraints for some $a_j$ (These are the energy constraints mentioned above.

14) $\frac{d(a_j)}{dt}=\frac{\frac{dE}{dt}-\int_0^{d} \sum_{i \neq j} c(y)\frac{d(a_i)}{dt} exp(-\lambda_i y)dy}{\int_0^{d}c(y)exp(-\lambda_j y)dy}$

15) $a_j=\frac{E-\int_0^{d} c(y)\sum_{i \neq j} a_i exp(-\lambda_i y)dy}{\int_0^{d} c(y) exp(-\lambda_i y)dy}$

Now the Euler Lagrange Equation is used to obtain the dynamics (differential equations) from the Lagrangian

16) $\left(\frac{\mathcal{L}(a_1,\ldots,a_n,\ldots,\dot a_1,\ldots,\dot a_n,\ldots;t)}{\partial a_i}\right) - \left(\frac{d}{dy} \frac{\mathcal{L}(a_1,\ldots,a_n, \ldots,;\dot a_1,\ldots,\dot a_n, \ldots;t)}{\partial \dot a_i}\right) = 0 \quad \textrm{for } i = 1, \dots, n, i \neq j$

This will give a second order differential equation and it might be necessary to use some linear algebra to rearrange the equation. You can convert this into a first order differential equation by using the Hamiltonian form. If only three basis functions are used then one gets “two mode model”. There are no restrictions on the number of bais functions used.

Like this:

Be the first to like this post.

September 1, 2009 - Posted by s243a | GCM (General Circulation Models)

7 Comments »

I saw a link for vertical ocean mixing posted on Lucia’s blog. It looked interesting:
http://oceanworld.tamu.edu/resources/ocng_textbook/chapter08/chapter08_05.htm

Also, I haven’t taken the above equations beyond this post but it looks like if I include my entropy condition in the Lagrangian then the equations become nonlinear but if I leave the entropy condition out the equations are linear.

Comment by s243a | September 1, 2009 | Reply
Just as a side note, I realized after posting this thread that the two box model is also an approximation to the heat equation. I discuss the case of overlapping boxes here:

Overlapping Boxes

Also of interest is a post I saw on the Crank-Nicolson method which is the standard technique for solving the heat equation:

Approximations used in Crank-Nicolson method for solving PDEs numerically

Comment by s243a | September 8, 2009 | Reply
I just wanted to share a good post with regards to modeling that I read over at the blackboard:

Bill Illis (Comment#19674) September 9th, 2009 at 10:51 pm

There has been a few points about volcanic forcing and the fast response. Effectively, volcanoes only have about one-third of the temperature impact as would be expected based on the decline in solar energy produced by the sulfate aerosols.

The climate models deal with this by building in an “efficacy of climate forcing” factor into all the different forcing. I guess one could think of this as a fudge factor since there should be no reason to expect a different response to each watt of change. A watt should be a watt. [I would rewrite the equations using the Stefan Boltzmann equations so that a watt is a watt but I guess that is for another day].

One can read a little more about this at this link but I can’t explain the rationale fully.

http://data.giss.nasa.gov/efficacy/

http://rankexploits.com/musings/2009/arthurs-case-2-i-think/

Comment by s243a | September 10, 2009 | Reply
I’m just going to post this link here for future reference. It is an old post I made about approximating nonlinear systems using the matrix exponential

http://www.physicsforums.com/showthread.php?t=325690

Comment by s243a | September 10, 2009 | Reply
Here is some good links on ocean dynamics I’ll want to look at later:
Ocean Dynamics

Comment by s243a | September 12, 2009 | Reply
[...] suggests my choice of a negative exponential basis in my post (Lagrangian Mechanics and the Heat Equation) was not two bad a choice although, Guasian functions will decay slightly faster then negative [...]

Pingback by Numeric Solutions to The Heat Equation « Earth Cubed | September 15, 2009 | Reply
I don’t know If I said it already but …I’m so glad I found this site…Keep up the good work I read a lot of blogs on a daily basis and for the most part, people lack substance but, I just wanted to make a quick comment to say GREAT blog. Thanks,

A definite great read..Tony Brown

Comment by Tony Brown | September 23, 2009 | Reply

Earth Cubed

Distributed Climate Science and Computing