2009 August « Earth Cubed

Laplace Transform of f(t) Related to smoothed f(t)?

When reading (Comment#18839) I started to wonder if there was a relationship between the Fourier Transform of a smoothed signal and the Laplace transform. I assumed there was a relationship (Comment#18854). After further derivation, I recommenced that if the goal is to derive the Laplace transform from the Fourier transom of the filtered signal:

1) The signal be properly windowed.

2) The FFT of the windowed Fourier Transform, needs to be compensated for the frequency effects that resulted from the low pass filter.

Weather it is a good idea to compute the Laplace transform from a windowed FFT of a filtered signal is outside of the scope of this thread (but feel free to comment bellow) .

The Laplace transform is given by:

1) $(\mathcal{L}f)(s) = \int_{0^-}^\infty e^{-st}f(t)\,dt$

The Fourier transform is given by:

2) $\hat{f}(\omega) := \int_{-\infty}^{\infty} f(x)\ e^{- \omega i t }\,dt,$

The Two Sided Laplace Transform is given by:

3) $\mathcal{B} \left\{f(t)\right\} = F(s) = \int_{-\infty}^{\infty} e^{-st} f(t) dt.$

Therefor the Fourier transform is the two sided Laplace transform evaluated at $s=i\omega$

Returning to the one sided Laplace transform:

4) $(\mathcal{L}f)(s) = \int_{a}^b e^{-st}f(t)\,dt= \int_{0^-}^\infty e^{-(\alpha+iw)t}f(t)\,dt$

5) $(\mathcal{L}f)(s) = \int_{a}^b e^{-st}f(t)\,dt= \int_{0^-}^\infty e^{-iw t}f(t)e^{- \alpha t}\,dt$

Integrating by Parts:

Let:

6) $u=e^{- iw t}$
7) $dv=f(t)e^{- \alpha t}\,dt$

8) $V=\int_{a}^tf(t')e^{-\alpha t'}\,dt'=e^{-\alpha t}\int_{t'=a}^tf(t')e^{\alpha (t-t')}\,dt'$

where the low pass filtered version of f(t):

9) $g(t)=f(t)*e^{\alpha t} =e^{-\alpha t} \int_{a}^tf(t')e^{\alpha (t-t')}\,dt'$

and is the convolution of f(t) and the impulse response of a filter (or atleast aproximatly so) with bandwidth $\alpha$ .

Plugging this result into integration by parts gives:

10) $(\mathcal{L}f)(s)=[e^{i \omega t}g(t)]_a^b-\int_a^b-i \omega e^{-i \omega t}e^{-\alpha t}g(t)dt$

or equivalently:

11) $(\mathcal{L}f)(s)=g(b)e^{-i\omega b}-g(a)e^{-i \omega a}+i\int_a^b \omega e^{-i \omega t}e^{-\alpha t}g(t)dt$

The first two terms show how the endpoints chosen effect the transfrom. These two terms will cancel for a given frequency if the distance between the endpoints is some multiple of the period. The last term is the Fouier transform of the smoothed function with the frequencies weighted by $\omega$ and using a windowing function $e^{-\alpha t}$
(note the multiple $i$ is there because the Fourier transform variable $\omega$ is the Laplace transform variable but rotated by 90 degrees.)

The effect of the windowing function is to smooth the frequency response. This is because multiplication in the time domain is equivalent to convolution in the frequency domain. The following Fourier transform relationship is useful (relationship 205):

12) $\mathcal{F}(e^{- \alpha t} u(t)) \,=\frac{1}{\sqrt{2 \pi} (\alpha + i \omega)}=\frac{1}{\sqrt{2 \pi (\alpha^2 + \omega^2)}}e^{(-i \ tan^{-1}(\alpha,\omega))}$

Note, that if a non causal filter was used for the smoothing the relationship is much simpler.

13) $\mathcal{F}({e}^{-a|t|}) \,=\sqrt{\frac{2}{\pi}} \cdot \frac{a}{a^2 + \omega^2}$

In both cases to properly deal with the end points the time shifting property of the Fourier transform is needed (relationship 102):

14) $f(x - t)\, =e^{- i a \omega} \hat{f}(\omega)\,$

Applying this property to the last two relationships gives:

15) $\mathcal{F}(e^{- \alpha (t-a)} (u(t-a)-u(t-b)) \,=\frac{e^{-ia \omega}-e^{\alpha a-ib\omega}}{\sqrt{2 \pi} (\alpha + i \omega)}$

16) $\mathcal{F}({e}^{-a|(t-b)|}) \,=\sqrt{\frac{2}{\pi}} \cdot \frac{a}{a^2 + \omega^2}e^{-i\omega\frac{1}{2}(a+b)}$

Strictly dealing with the case where a causal filter is used and applying the rule for the Fourier transform of a convolution (Rule 109) we obtain:

17) $(\mathcal{L}f)(s)=g(b)e^{-i\omega b}-g(a)e^{-i \omega a}+i\omega\frac{1}{\sqrt{2\pi}}\mathcal{F}(g(t))*\left( \frac{e^{-ia \omega}-e^{\alpha a-ib\omega}}{\sqrt{2 \pi} (\alpha + i \omega)} \right)$

of equivalently:

18) $\mathcal{L}(f(y))=g(b)e^{-i\omega b}-g(a)e^{-i \omega a}+i\omega\frac{1}{\sqrt{2\pi}}\mathcal{F}(g(t))*\left( \frac{e^{-ia \omega}-e^{\alpha a-ib\omega}}{\sqrt{2 \pi(\alpha^2 +\omega^2)}}e^{ -i tan^{-1}(\alpha,\omega) } \right)$

Some Comments:

1) If $\alpha$ is negative the system is causal, and the filtered version $g(t)$ of the signal $f(t)$ will be causal.

2) Computing the smoothed signal does not save any computations with regards to computing the Laplace transform.

3) The derivation seems to show that their is a relationship between Laplace trancform and a windowed Fouier transform of the filtered signal.

4) To compute the Laplace transform based on the orginal signal use equation (5). To compute it based on the filtered signal use equation (11).

August 30, 2009 Posted by s243a | Math | 8 Comments

Coriolis Forces

A derivation for coriolis forces can be found be found on the wikipedia page for fictitious forces.

In general for an accelerating reference frame in rectangular coordinates the factious forces are given by:

$\mathbf{F}_{\mbox{fictitious}} =-m\ \mathbf{a}_{AB} -2m\ \sum_{j=1}^3 v_j \ \frac{d \mathbf{u}_j}{dt} - m\ \sum_{j=1}^3 x_j \ \frac{d^2 \mathbf{u}_j}{dt^2}\ .$

Where:

he first term is the Coriolis force, the second term is the centrifugal force, and the third term is the Euler force. When the rate of rotation doesn’t change, as is typically the case for a planet, the Euler force is zero.

Looking specifically at the Coriolis force:

$- 2 m \boldsymbol \omega \times \mathbf{v}$

which gives in (east-west,north-sourth, height) coordinates:

$\boldsymbol{ \Omega} = \omega \begin{pmatrix} 0 \\ \cos \varphi \\ \sin \varphi \end{pmatrix}\ ,$ $\boldsymbol{ v} = \begin{pmatrix} v_e \\ v_n \\ v_u \end{pmatrix}\ ,$

$\boldsymbol{ a}_C =-2\boldsymbol{\Omega \times v}= 2\,\omega\, \begin{pmatrix} v_n \sin \varphi-v_u \cos \varphi \\ -v_e \sin \varphi \\ v_e \cos\varphi\end{pmatrix}\ .$

Where $\phi$ is the latitudinal coordinate (equator=zero latitude).

In general the cross product for a coordinate with orthonormal direction vectors is given by:

$\boldsymbol{\Omega \times v} = \begin{vmatrix} \boldsymbol{i}&\boldsymbol{j}&\boldsymbol{k} \\ \Omega_x & \Omega_y & \Omega_z \\ v_x & v_y & v_z \end{vmatrix}\ = \begin{pmatrix} \Omega_y v_z - \Omega_z v_y \\ \Omega_z v_x - \Omega_x v_z \\ \Omega_x v_y - \Omega_y v_x \end{pmatrix}\ ,$
http://en.wikipedia.org/wiki/Coriolis_effect#Formula

since the basis direction vectors are orthogonal in hopkins and simmons coordinates write write:

$\boldsymbol{\Omega \times v} = \begin{vmatrix} \boldsymbol{e_{\sigma}}&\boldsymbol{e_{\mu}}&\boldsymbol{e_{\lambda}} \\ \Omega_{\sigma} & \Omega_{\mu} & \Omega_{\lambda} \\ W & V & U \end{vmatrix}\ = \begin{pmatrix} \Omega_{\mu} U - \Omega_{\lambda} V \\ \Omega_{\lambda} W - \Omega_{\sigma} U \\ \Omega_{\sigma} V - \Omega_{\mu} W \end{pmatrix}\ ,$

(note with regards to weather the system is right handed we can choose the direction of the logitudanal cordinate $\lambda$ to make it right handed.)

Just to recall from the post (Hoskins and Simmons (1974) Coordinate System):

$\mu = sin( \theta )$ where theta is the latitude.
$\sigma = pressure/P_*$ Where $P_*$ is the surface pressure and $\sigma$ is the vertical coordinate.
$\lambda$ is the longitude.

Additionally:

U is the longitudinal component of the velocity
V is the latitudinal component of the velocity
W is the vertical component of the veolicty (not used in Hopkins and Simmons 1974)

Now the angular velocity of the earth in Hopkins and Simmons is given by:

$\boldsymbol{ \Omega} = \omega \begin{pmatrix} \cos \varphi \\ \sin \varphi \\ 0 \end{pmatrix}\ = \omega \begin{pmatrix} \mp \sqrt{1-\mu^2} \\ \mu \\ 0 \end{pmatrix}\,$ $\boldsymbol{ v} = \begin{pmatrix} W \\ V \\ U \end{pmatrix}\ ,$

Where the sign of $\pm$ is positive for the northern hemisphere and negative for the southern hemisphere.

Therefore:

$\boldsymbol{\Omega \times v} = \begin{pmatrix} \mu U \\ \pm \sqrt{1-\mu^2} U \\ \mp \sqrt{1-\mu^2} V - \mu W \end{pmatrix} \$

Some comments:

The result obtained is essentialy the same result that one would get if, they took the(east west,north south, altitude) coordinate system and replaced $\phi$ with $asin(\mu)$ .

The only differences are the order and sign of the components. These are the only differences because both coordinate contain the same unit vectors. In my example of a Hopkins and Simon’s like coordinate system I used a different order for the components then was used in my example for the (east-west, north south altitude) coordinate system. This will effect the sign in the cross product.

I wrote the z component of the angular velocity as $\mp \sqrt{1-\mu^2}$ to emphasize that the positive direction for the z component in Simpons coordinate system is downward. However, the actual angular rotation of the earth in simons coordinate system still have a postive $\mu$ component depending on the which direction is defined as positive for the longitudinal coordinate.

The order which we specify the coordinates determines the right handedness of the coordinate system. Therefore, righthandedness is not inherently a geometric property because it depends on the order of the coordinates. For instance, in standard Cartesian coordinates $e_z \times e_y = -e_x$

In our case the first coordinat,e $\sigma$ , was specified in the downward direction, our second coordinate, $\mu$ , points south, now using the right hand rule means that gives the positive direction for the third coordinate $\lambda$ in the east direction.

It is for these reasons that differences can arrise, and therefore it is very important when doing cross products to clearly express the postive direction of the coordinate unit vectors and the order of the coordinates.

August 29, 2009 Posted by s243a | GCM (General Circulation Models) | 2 Comments

Vector Operations in Hoskins and Simmons Coordinates

In my post Hoskins and Simmons (1974) Coordinate System, I derived the following scaling quantities which will be used to derive the vector operations of Grad, Div and Curl in Hoskins Coordinate system.

$h_{\mu}= \left| \frac{r}{\sqrt{1-\mu^2}} \right|$
$h_{\lambda}=r \mu$
$h_{\sigma} =-\sigma \frac{mg}{RT}$

The coordinates in Hoskins coordinate system are dimensionless . (see nondimentionalization of Navier Stokes).

The gradient is defined as (see lectures on coordinate transforms):

$\nabla \phi = \frac{\vec e_{\mu}}{h_{\mu}} \frac{\partial \phi}{\partial \mu} + \frac{\vec e_{\lambda}}{h_{\lambda}}\frac{\partial \phi}{\partial \lambda}+\frac{\vec e_{\sigma}}{h_{\sigma}}\frac{\partial \phi}{\partial \sigma}=\left(\vec e_{\mu}\right)\frac{1-\mu^2}{r} \frac{\partial \phi}{\partial \mu} + \left(\vec e_{\lambda}\right)\frac{1}{r \mu}\frac{\partial \phi}{\partial \lambda}-\left(\vec e_{\sigma}\right)\frac{RT}{\sigma mg}\frac{\partial \phi}{\partial \sigma}$

The divergence is defined as:

$\nabla \cdot = \frac{1}{h_{\mu}} \frac{\partial }{\partial \mu} + \frac{1}{h_{\lambda}}\frac{\partial}{\partial \lambda}+\frac{1}{h_{\sigma}}\frac{\partial}{\partial \sigma}=\frac{1-\mu^2}{r} \frac{\partial}{\partial \mu} + \frac{1}{r \mu}\frac{\partial}{\partial \lambda}-\frac{RT}{\sigma mg}\frac{\partial}{\partial \sigma}$

The curl is defined by:

$\nabla \times \vec A = \frac{1}{h_{\sigma}h_{\mu}h_{\lambda}} \left| \begin{array}{ccc} \vec e_{\sigma} & \vec e_{\mu} & \vec e_{\lambda} \\ \frac{\partial}{\partial \sigma} & \frac{\partial}{\partial \mu} & \frac{\partial}{\partial \lambda} \\ A_{\sigma}h_{\sigma} & A_{\mu}h_{\mu} & A_{\lambda}h_{\lambda} \end{array} \right|$

(note, the direction of the longitudinal coordinate $\lambda$ is defined to obey the right hand rule)

This gives for the components

$\left[ \nabla \times \vec A \right]_{\sigma}=\frac{RT}{r \mu \sigma mg}\left(\sigma \frac{mg}{RT}\frac{\partial A_{\lambda}}{\partial \mu}+ r \mu\frac{\partial A_{\mu}}{\partial \lambda} \right)$

$\left[ \nabla \times \vec A \right]_{\mu}= -\frac{\sqrt{1-\mu^2}RT}{r\sigma mg}\left( \frac{r}{\sqrt{1-\mu^2}}\frac{\partial A_{\sigma}}{\partial \lambda}+\sigma \frac{mg}{RT} \frac{\partial A_{\lambda}}{\partial \sigma} \right)$

$\left[ \nabla \times \vec A \right]_{\lambda}= \frac{\sqrt{1-\mu^2}RT}{r^2 \mu }\left( r \mu\frac{\partial A_{\mu}}{\partial \sigma}-\frac{r}{\sqrt{1-\mu^2}} \frac{\partial A_{\sigma}}{\partial \mu} \right)$

Which Simplifies to:

$\left[ \nabla \times \vec A \right]_{\sigma}= \frac{1}{r\mu}\frac{\partial A_{\lambda}}{\partial \mu}+ \frac{RT}{\sigma m g}\frac{\partial A_{\mu}}{\partial \lambda}$

$\left[ \nabla \times \vec A \right]_{\mu}= -\frac{RT}{\sigma m g}\frac{\partial A_{\sigma}}{\partial \lambda}- \frac{\sqrt{1-\mu^2}}{r} \frac{\partial A_{\lambda}}{\partial \sigma}$

$\left[ \nabla \times \vec A \right]_{\lambda}= \frac{\sqrt{1-\mu^2}RT}{r}\frac{\partial A_{\mu}}{\partial \sigma}-\frac{RT}{r\mu} \frac{\partial A_{\sigma}}{\partial \mu}$

August 29, 2009 Posted by s243a | Uncategorized | 1 Comment

API/Object Viewers/Memory Mapping/

The more code a programmer can reuse the more efficient they can be. In windows this could mean reusing com/ole components and other APIs. Here are two useful programs for viewing APIs:

OLE/COM Object Explorer 1.1
http://www.softpedia.com/progDownload/OLE-COM-Object-Explorer-Download-42531.html

Windows API Viewer
http://www.activevb.de/rubriken/apiviewer/index-apiviewer.html

I was inquiring about how to manage the transfer of large amounts of data between programs and I was pointed to two interesting concepts:

All modern operating systems include a facility called “memory mapping,” which maps a range of addresses in the program’s virtual address space to a file. If you read from those addresses, you’ll get data from the file. It is up to the operating system to determine whether to load the data into RAM all at once, or to read it from the disk in chunks as necessary.
…….
If you’re trying to share large amounts of memory between two programs running on the same computer, you should note that all modern operating systems provide mechanisms for shared memory. These shared memory segments can be mapped into the virtual address space of multiple programs simultaneously. Two or more programs can read or write to the shared memory exactly as if it were normal, private memory. (But you should include some thread-safety mechanisms, like mutexes, to make sure your programs won’t step on each other’s toes.)

If you’re trying to share large amounts of memory between programs running on separate computers, use MPI or some other multi-processing library.

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

Here is what wikipedia has to say about memory maps:

The primary benefit of memory mapping a file is increased I/O performance, especially when used on small files. Accessing memory mapped files is faster than using direct read and write operations for two reasons. Firstly, a system call is orders of magnitude slower than a simple change of program’s local memory. Secondly, in most operating systems the memory region mapped actually is the kernel’s file cache, meaning that no copies need to be created in user space. Using system calls would inevitably involve the time consuming operation of memory copying.

Certain application level memory-mapped file operations also perform better than their physical file counterparts. Applications can access and update data in the file directly and in-place, as opposed to seeking from the start of the file or rewriting the entire edited contents to a temporary location. Since the memory-mapped file is handled internally in pages, linear file access (as seen, for example, in flat file data storage or configuration files) requires disk access only when a new page boundary is crossed, and can write larger sections of the file to disk in a single operation.

A possible benefit of memory-mapped files is a “lazy loading”, thus using small amounts of RAM even for a very large file. Trying to load the entire contents of a file that is significantly larger than the amount of memory available can cause severe thrashing as the operating system reads from disk into memory and simultaneously pages from memory back to disk. Memory-mapping may not only bypass the page file completely, but the system only needs to load the smaller page-sized sections as data is being edited, similarly to demand paging scheme used for programs.

http://en.wikipedia.org/wiki/Memory-mapped_file#Benefits

The windows utility to do this is called CreateFileMapping
http://msdn2.microsoft.com/en-us/library/aa366537.aspx

As for multi-processing library the following was recommended:
http://scv.bu.edu/documentation/tutorials/MPI/MPI_text.html

I haven’t found much but the following wikipedia link seems relevant:
http://en.wikipedia.org/wiki/Cluster_(computing)

August 29, 2009 Posted by s243a | Computer Science and Modeling | 2 Comments

Defining a Microsoft access Datasource

The more I think about implementation of a GCM the more I think about trying to make the components as independent as possible. I like databases for storing information because they give an orderly framework that is easy to access but they can be slow.

A database is simply a binary file. If we know how information is stored in the file we can optimize ways of retrieving it that could be much faster then a standard database. However, it is nice to have tools that are easy to use for sorting though and organizing information like Microsoft access.

Remember that tables in microsoft access do not need to be stored in the microsoft access database. For instance, you can create linked tables. These linked tables can come for instance from other databases or spreadsheets. All that is really necessary is for Microsoft access to know how to read this information as a table.

This is where I believe data access objects come in. I’ll have to look at them further but I’m hoping that they can be used as a way to use an arbitrary source of information as a table. The following linking gives information on how to create Data Access objects in Microsoft Access.

http://msdn.microsoft.com/en-us/library/ms243192.aspx

I’ll write more once I read about this.

August 28, 2009 Posted by s243a | Computer Science and Modeling | Leave a Comment

Fractal Modeling of Turbulence

For now this is just a place holder to discuss topics about fractals. For related topics see my post on Kolmogorov’s Turbulance. I would like to present the quote though to illustrate the difficulty of using numerical methods to solve Naiver Stokes equations:

1.1. Statement of the problem Many flows of interest in science and engineering display complex spatial and temporal structures (eddies) spanning a wide range of scales. The ratio between the largest (L) and smallest ( $\eta$ ) scale can easily exceed $10^4$ in typical engineering applications, and can be as high as $10^6$ or higher in geophysical applications. Since the nonlinear interaction between eddies of different sizes eludes even the most sophisticated analytical approaches, one must resort to either extensive experimentation or direct numerical simulation (DNS) of the governing equations. The latter approach has gained strength by the rapid increase in the power of digital computers during the past 20 years. Despite this fact, DNS of flows for which the ratio $latex L/ \eta $ is much larger than $10^2$ are still prohibitive

http://www.me.jhu.edu/meneveau/pdf-papers/ScottiMeneveau99.pdf

More papers on fractals and turbulance can be found here:

http://www.me.jhu.edu/~meneveau/pubs-fractals.html

August 27, 2009 Posted by s243a | Fractals, Turbulence | 3 Comments

Hoskins and Simmons (1974) Coordinate System

sHoskins and Simmons (1974) Coordinate System

Hoskin and Simon used the following coordinate system:

$\mu = sin( \theta )$ where theta is the latitude.

$\sigma = pressure/P_*$ Where $P_*$ is the surface pressure and $\sigma$ is the verticle coordinate.

$\lambda$ is the longitude.

This coordinate system is also used in the following model which was recomended to me as a simple model:

As the complexity of general circulation models has been and still is growing considerably, it is not surprising that, for both education and research, models simpler than those comprehensive GCMs at the cutting edge of the development, are becoming more and more attractive. These medium complexity models do not simply enhance the climate model hierarchy. They support understanding atmospheric or climate phenomena by simplifying the system gradually to reveal
the key mechanisms. They also provide an ideal tool kit for students to be educated and to teach themselves, gaining practice in model building or modeling. Our aim is to provide such a model of intermediate complexity for the university environment: the PlanetSimulator. It
can be used for training the next GCM developers, to support scientists to understand climate processes, and to do fundamental research.

http://www.mi.uni-hamburg.de/216.0.html?&L=1

http://www.mi.uni-hamburg.de/fileadmin/files/forschung/theomet/planet_simulator/downloads/PS_ReferenceGuide.pdf

In order to use the coordinates given by Hoskins and Simmons in the Navier Stokes equations it is necissary to find the operations for Div Grad and Curl in this coordinate system. You can get a free PDF here:

http://www.scribd.com/doc/2590597/Lectures-on-Transformation-of-Coordinates

Which shows how to derive these operations for any coordinate system. These derivations express the apove operations in terms of scaling quantities. The scaling quantitys are the magnitude of the derivative of the change in postion, with respect to a change in one of the coordinates in the new coordinate system.

In our case for the coordinate $\mu$ the scalling quanity is:

$h_{\mu}=\left| \frac{ \partial \vec r}{\partial \mu} \right| = \left| \frac{ \partial (r \theta)}{\partial \mu} \right|$

because when you change latitude the distance moved is proportional to the radius multiplied by the change in latitudinal (measured in radians).

Now expressing the latitude in terms of $\mu$
$h_{\mu}=\left| \frac{ \partial (r \ sin^{-1}( \mu ))}{\partial \mu} \right| = \frac{r}{\sqrt{1-\mu^2}}$

Now for the $\lambda$ which is the longitude coordinate, $h_{\lambda}$ depends on how you define zero latitude. The distance changed for a given change in longitude is the product of the change in logitude (measured in radians) multiplied by the distance to the earths axis of rotation. If zero lattitude is the center of the earth then this distance is $r cos(\theta)$ if zero latitude is the north pole (typical in spherical coordinates) then this distance is $r sin(\theta)$ .

We will see that taking zero latitude as the north pole will be simpler.

$h_{\lambda}=\left| \frac{ \partial \vec r}{\partial \lambda} \right| = \left| \frac{ \partial (r \lambda sin(\theta))}{\partial \lambda} \right|=r \ sin(\theta)$

But $\theta=sin^{-1}(\mu)$ Therefore:

$h_{\lambda}=r \mu$

The vertical coordinate $\sigma$ is the ratio of the hydrostatic pressure to the hydrostatic pressure at sea level.

In Hoskins & Simmons the hydrostatic pressure used in the coordinate ignores the lapse rate. In some models variations on this are used. For instance in a model used by Judith Curry this coordinate is somewhat adjusted based on the surface topology of the earth. (see my post: hybird Coordinate System)

hybird coordinate system

In Curries model this is done to try and keep wind flow on the horizontal plane. Anyway, the equation for hydrostatic equilibrium is:

$\frac{dp}{dh} = - \frac {mpg}{RT}$

Now:

$\frac{\partial \sigma}{\partial h}=\frac{\partial \sigma}{\partial P}\frac{\partial P}{\partial h}=\frac{1}{P_*}\frac{\partial P}{h}=-\frac{1}{P_*} \frac {mpg}{RT} =-\sigma \frac{mg}{RT}$

Hmmmmmmm: From what I’ve done the only way I can see to make the constants $\frac{mg}{RT}$ disappear is to set $T=\frac{mg}{R}$

Here is what Hoskins and Simmons have to say about temperature scaling and such.

all non-dimensionalized using as length scale the radius of the planet, a, ; as time scale the reciprocal of its angular velocity $\Omega$ ; temperature scale $a^2\Omega^2/R$ (R as the gas constant) and pressure scale $P_o=1$ Bar.

In http://www.mi.uni-hamburg.de/fileadmin/files/forschung/theomet/planet_simulator/downloads/PS_ReferenceGuide.pdf $T_o$ was taken to be 250K. Some numbers of rellevence:

Gas constant of air: 287.058 J kg−1 K−1

and for the molar mass of air:

Standard temperature =273.16k

Source for the molar mass of air:

Standard pressure=1 Atm=1013.250*10^2 Pa.
The molar mass of air at standard temperature and pressure is 28.964*10^-3 kg/mole. Molar

http://www.blurtit.com/q139778.html

Standard Gravity: 9.80665 m/s2

August 27, 2009 Posted by s243a | GCM (General Circulation Models), Math | 3 Comments

Deriving The Vorticy Equation

As I mentioned here wikipedia has two good links on the vorticity equation:

http://en.wikipedia.org/wiki/Vorticity_equation
http://en.wikipedia.org/wiki/Barotropic_vorticity_equation

At the moment the derivation seems to be missing but there is enough information provided to do it:

From wikipedia’s link Derivation of Navier Stokes Equations the following equation is given:

$\rho\frac{D\mathbf{v}}{D t} = -\nabla p + \nabla \cdot\mathbb{T} + \mathbf{f}$

Dividing each side by the density you get roughly form of stokes equations given here (Note, we have $\vec B$ instead of $\rho \Vec B$ . We get this by letting B be the acceleration per unit mass, so the vorce times unit mass is the acceleration multiplied by the density. Then when you divide both thids of the above by the density we are left with just $\vec B$ this will work out to the correct result in the end.

$\frac{D \vec V}{D t} = \frac{\partial \vec V}{\partial t} + (\vec V \cdot \vec \nabla) \vec V = - \frac{1}{\rho} \vec \nabla p + \vec B + \frac{\vec \nabla \cdot \underline{\underline{\tau}}}{\rho}$

Now the vorticy equation is derived by taking the cross product of both sides. This will give on the Left hand side of the equation:

$LHS=\nabla \times \frac{D \vec V}{D t} = \nabla \times \left( \frac{\partial \vec V}{\partial t} + (\vec V \cdot \vec \nabla) \vec V \right)$

$LHS= \frac{\partial \vec \omega}{\partial t} + \nabla \times ( \vec V \cdot \vec \nabla ) \vec V$

using the identity (see the thread: Vector Identities and Derivations)

$\vec V \cdot \vec \nabla \vec V = \vec \nabla (\tfrac{1}{2} \vec V \cdot \vec V) - \vec V \times \vec \omega$

We get:

$LHS= \frac{\partial \vec \omega}{\partial t} + \nabla \times \left( \vec \nabla (\tfrac{1}{2} \vec V \cdot \vec V) - \vec V \times \vec \omega \right)$

since the curl of the divergence is equal to zero this reduces to:

$LHS= \frac{\partial \vec \omega}{\partial t} - \nabla \times \left( \vec V \times \vec \omega \right)$

using the identity:

$\vec \nabla \times (\vec V \times \vec \omega ) = -\vec \omega (\vec \nabla \cdot \vec V) + (\vec \omega \cdot \vec \nabla ) \vec V - (\vec V \cdot \vec \nabla) \vec \omega$

We arrive at:

$LHS= \frac{\partial \vec \omega}{\partial t} + \vec \omega (\vec \nabla \cdot \vec V) - (\vec \omega \cdot \vec \nabla ) \vec V + (\vec V \cdot \vec \nabla) \vec \omega$

Now combining this with the cross product of the right hand side one obtains:

$\frac{\partial \vec \omega}{\partial t} + \vec \omega (\vec \nabla \cdot \vec V) - (\vec \omega \cdot \vec \nabla ) \vec V + (\vec V \cdot \vec \nabla) \vec \omega = - \nabla \times \left( \frac{1}{\rho} \vec \nabla p \right) + \nabla \times \vec B + \nabla \times \frac{\vec \nabla \cdot \underline{\underline{\tau}}}{\rho}$

rearanging:

$\frac{\partial \vec \omega}{\partial t} + (\vec V \cdot \vec \nabla) \vec \omega = (\vec \omega \cdot \vec \nabla ) \vec V - \vec \omega (\vec \nabla \cdot \vec V) - \nabla \times \left( \frac{1}{\rho} \vec \nabla p \right) + \nabla \times \vec B + \nabla \times \frac{\vec \nabla \cdot \underline{\underline{\tau}}}{\rho}$

This is nearly equivalent to the form of the vorticity equation shown in Wikipedia except for this term:

$- \nabla \times \frac{1}{\rho} \vec \nabla p$

The following identity is needed:

$\nabla \times (\psi\mathbf{A}) = \psi\nabla \times \mathbf{A} - \mathbf{A} \times \nabla\psi$

Therefore:

$- \nabla \times \frac{1}{\rho} \vec \nabla p = \frac{1}{\rho} \nabla \times \vec \nabla p - \vec \nabla p \times \nabla \frac{1}{\rho}$

but since the curl of a gradient is equal to zero:

$- \nabla \times \frac{1}{\rho} \vec \nabla p = - \vec \nabla p \times \nabla \frac{1}{\rho}$

Now applying the chain rule:

$- \nabla \times \frac{1}{\rho} \vec \nabla p = - \vec \nabla p \times \frac{1}{\rho^2} \nabla \rho$

Reversing the order of the cross product changes the sign. Consequently:

$- \nabla \times \frac{1}{\rho} \vec \nabla p = \frac{1}{\rho^2} \nabla \rho \times \vec \nabla p$

substituting this result back into the vorticity equation gives:

$\frac{\partial \vec \omega}{\partial t} + (\vec V \cdot \vec \nabla) \vec \omega = (\vec \omega \cdot \vec \nabla ) \vec V - \vec \omega (\vec \nabla \cdot \vec V) + \frac{1}{\rho^2} \nabla \rho \times \vec \nabla p + \nabla \times \vec B + \nabla \times \frac{\vec \nabla \cdot \underline{\underline{\tau}}}{\rho}$

The following simplifications will be useful depending on the application:

# In case of conservative force, $\vec \nabla \times \vec B = 0$
# For barotropic fluid, $\vec \nabla \rho \times \vec \nabla p = 0$ . This is also true for a constant density fluid where $\vec \nabla \rho = 0$

# For inviscid fluids, $\underline{\underline{\tau}} = 0$ .

Vector Identities and Derivations

August 25, 2009 Posted by s243a | Navier Stokes | 3 Comments

Kolmogorov’s Theory of Turbulence

It is impractical in climate models to give sufficient resolution in order to capture all of the inertial dynamics in the system. Inertial dynamics will dominate at large scales while on small scales viscous forces will despite the energy associated with these dynamics.

The Reynolds number gives the ratio of inertial forces to viscous forces. When viscous forces dominate it is known as laminar flow. Laminar flow occurs when Reynolds number is less then 10. If velocities is held constant and the length scale decreases then so does reynolds number and for this small region the flow will approach Laminar flow and average shear forces will be dominated by the viscosity.

When viscous forces dominate the shear forces are dominated by the momentum exchange though partial diffusion, when viscous forces don’t dominate then the bulk flow of fluid can give a momentum exchange that acts like an effective increase in viscosity. Therefore when making large scale approximations it is important to know how the sub scale dynamics effect large scale properties such as the average shear forces.

One way we might do this is though Kolmogrov’s Theory of Turbulence:

A turbulent flow is characterized by a hierarchy of scales through which the energy cascade takes place. Dissipation of kinetic energy takes place at scales of the order of Kolmogorov length η, while the input of energy into the cascade comes from the decay of the large scales, of order L. These two scales at the extremes of the cascade can differ by several orders of magnitude at high Reynolds numbers. In between there is a range of scales (each one with its own characteristic length r) that has formed at the expense of the energy of the large ones. These scales are very large compared with the Kolmogorov length, but still very small compared with the large scale of the flow (i.e. $\eta \ll r \ll L$ ). Since eddies in this range are much larger than the dissipative eddies that exist at Kolmogorov scales, kinetic energy is essentially not dissipated in this range, and it is merely transferred to smaller scales until viscous effects become important as the order of the Kolmogorov scale is approached. Within this range inertial effects are still much larger than viscous effects, and it is possible to assume that viscosity does not play a role in their internal dynamics (for this reason this range is called “inertial range”).

http://en.wikipedia.org/wiki/Turbulence#Kolmogorov.27s_Theory_of_1941

Some other relevant links:

http://www.sjsu.edu/faculty/watkins/kolmo.htm
http://en.wikipedia.org/wiki/Kolmogorov_microscales

http://scienceworld.wolfram.com/physics/Turbulence.html

On a related curiosity, I am wondering how knowing the flow conditions at a discrete number of points effects the range of possible solutions to Navier Stokes Equations. I’ve posted a bit on this in the next two links:

http://www.climateaudit.org/phpBB3/viewtopic.php?f=4&t=760
http://www.physicsforums.com/showthread.php?t=331713

August 24, 2009 Posted by s243a | Turbulence | 14 Comments

Vector Identities and Derivations

To derive the prognostic equations in a GCM usually requires the use of vector identities. Here is a list of some links which give vector identities:

http://en.wikipedia.org/wiki/List_of_vector_identities
http://en.wikipedia.org/wiki/Product_rule
http://en.wikipedia.org/wiki/Vector_calculus_identities
http://mathworld.wolfram.com/VectorDerivative.html
http://planetphysics.org/encyclopedia/VectorIdentities.html
http://geo.phys.spbu.ru/~runov/SPIntro/SPIntro_15.pdf

I’ve been searching for the derivations of these identities and have found an excellent source:

Vector algebra is a powerful and needful tool for Physics but unfortunately, due to lack of mathematical skills, it becomes misleading for first undergraduate courses of science and engineering studies. Standard vector identities are usually proved using Cartesian components or geometrical arguments, accordingly. Instead, this work presents a new teaching strategy in order to derive symbolically vector identities without analytical expansions in components, either explicitly or using indicial notation. This strategy is mainly based on the correspondence between three-dimensional vectors and skew-symmetric second-rank tensors. Hence, the derivations are performed from skew tensors and dyadic products, rather than cross products. Some examples of skew-symmetric tensors in Physics are illustrated.

http://www.citeulike.org/user/pak/article/4524046
http://arxiv.org/PS_cache/arxiv/pdf/0904/0904.1814v1.pdf

This will require some understanding of some concepts of tensor algebra. Here are some helpfull links of concepts that will be needed:

http://www.foamcfd.org/Nabla/guides/ProgrammersGuidese3.html
http://en.wikipedia.org/wiki/Levi-Civita_symbol
http://en.wikipedia.org/wiki/Dyadics
http://en.wikipedia.org/wiki/Pseudovector

As a side note, the prognostic equations are typically derived from the vorticy equation. See:

http://en.wikipedia.org/wiki/Vorticity_equation
http://en.wikipedia.org/wiki/Barotropic_vorticity_equation

August 23, 2009 Posted by s243a | Math | 2 Comments

« Previous Entries