The neutral density ( ) or empirical neutral density is a density variable used in oceanography, introduced in 1997 by David R. Jackett and Trevor McDougall. It is function of the three state variables (salinity, temperature, and pressure) and the geographical location (longitude and latitude) and it has the typical units of density (M/V). The level surfaces of form the “neutral density surfaces”. These are widely regarded as the most natural layer interfaces stratifying the deep ocean circulation, along which the strong lateral mixing in the ocean occurs, although this has yet to be rigorously established. These surfaces are widely used in water masses analysis. Neutral density is a density variable that depends on the particular state of the ocean, and hence a function of time as well, although this is often ignored. In practice, its construction for a given hydrographic observation is achieved by means of a computational code (available for Matlab and Fortran), that contains the computational algorithm developed by Jackett and McDougall. Use of this code is currently restricted to the present day ocean.
McDougall and Jackett  demonstrated that the normal to the neutral surfaces is in the direction of , where S is the salinity, is the potential temperature, the thermal expansion coefficient and the saline concentration coefficient. Thus, neutral surfaces are defined as the surfaces everywhere perpendicular to the vector . For such a surface to exist, its helicity H must be zero; if this condition is respected, a scalar exists and it is the one which satisfies the following formula:
where b is an integrating scalar factor, which is function of space, and a residual term that is usually different from zero.
This formula represents a coupled system of first-order partial differential equations, that has to be solved to obtain the desired value of . The solutions of (1 ) can be obtained by using numerical techniques.
In the real ocean, the condition of helicity equal to zero is not generally satisfied exactly. Therefore, and because of the non-linear terms in the equation of state, it is impossible to create analytically a Well-defined neutral density surface. There will always be flow through the calculated surfaces, because of the presence of a neutral helicity.
Therefore, it is possible to obtain only a best-fit approximate neutral surface, through which there is no flow of major proportions and along which it is generally accepted that flow takes place. is a Well-defined function and Jackett and McDougall demonstrated that the inaccuracy due to the not exact neutrality is below the present instrumentation error in density. Neutral density surfaces stay within a few tens meters of an ideal surface anywhere in the world.
For how has been defined, neutral density surfaces can be considered the continuous analog of the commonly used potential density surfaces, which are defined over various discrete values of pressures (see for example  and ).
Given the spatial dependence of the neutral density, its calculation requires the knowledge of the spatial distribution of temperature and salinity in the ocean. Therefore, the definition of has to be linked with a global hydrographic dataset, based on the climatology of the world’s ocean (see World Ocean Atlas and ). In this way, the solution of (1 ) provides values of for a referenced global dataset. The solution of the system for a high resolution dataset would be computationally very expensive. In this case, the original dataset can be sub-sampled and (1 ) can be solved over a more limited set of data.
Algorithm for the computation of neutral surfaces using 
Jackett and McDougall constructed the variable using the data in the “Levitus dataset”. As this dataset consist of measurements of S and T at 33 standard depth levels at a 1° resolution, the solution of (1 ) for such a large dataset would be computationally very expensive. Therefore, they sub-sampled the data of the original dataset onto a 4°x4° grid and solved (1 ) on the nodes of this grid. The authors suggested to solve this system by using a combination of the method of characteristics in nearly 85% of the ocean (the characteristic surfaces of (1 ) are neutral surfaces along which is constant) and the finite differences method in the remaining 15%. The output of these calculations is a global dataset labeled with values of . The field of values resulting from the solution of the differential system (1 ) satisfies (1 ) an order of magnitude better (on average) than the present instrumentation error in density.
The labeled dataset is then used to assign values to any arbitrary hydrographic data at new locations, where values are measured as a function of depth by interpolation to the four closest points in the Levitus atlas.
Practical computation of 
The Neutral Density code comes as a package of Matlab or as a Fortran routine. It enables the user to fit neutral density surfaces to arbitrary hydrographic data and just 2 MBytes of storage are required to obtain an accurately pre-labelled world ocean.
Then, the code permits to interpolate the labeled data in terms of spatial location and hydrography. By taking a weighted average of the four closest casts from the labeled data set, it enables to assign values to any arbitrary hydrographic data.
Advantages of using the neutral density variable
Comparisons between the approximated neutral surfaces obtained by using the variable and the previous commonly used methods to obtain discretely referenced neutral surfaces (see for example Reid (1994), that proposed to approximate neutral surfaces by a linked sequence of potential density surfaces referred to a discrete set of reference pressures) have shown an improvement of accuracy (by a factor of about 5)  and an easier and computationally less expensive algorithm to form neutral surfaces. A neutral surface defined using differs only slightly from an ideal neutral surface. In fact, if a parcel moves around a gyre on the neutral surface and returns to its starting location, its depth at the end will differ by around 10m from the depth at the start. If potential density surfaces are used, the difference can be hundreds of meters, a far larger error.
- Jackett, David R., Trevor J. McDougall, 1997: A Neutral Density Variable for the World's Oceans. J. Phys. Oceanogr., 27, 237–263
- McDougall, T. J. and D. R. Jackett, 1988: On the helical nature of neutral surfaces. Progress in Oceanography, Vol. 20, Pergamon, 153–183
- Jackett, David R., Trevor J. McDougall, 1997: A Neutral Density Variable for the World's Oceans. J. Phys. Oceanogr., 27, page 238
- Klocker et all., 2007, “Diapycnal motion due to neutral helicity”
- Jackett, David R., Trevor J. McDougall, 1997: A Neutral Density Variable for the World's Oceans. J. Phys. Oceanogr., 27, page 239
- "Introduction to Physical Oceanography : Chapter 6 - Temperature, Salinity, and Density - Density, Potential Temperature, and Neutral Density". Oceanworld.tamu.edu. Retrieved 2012-11-16.
- Montgomery, R. B., 1938: Circulation in the upper layers of the southern North Atlantic, Pap. Phys. Oceanogr. Meteor., 6(2), 55 pp.
- Reid, J. L., 1994: On the total geostrophic circulation of the North Atlantic Ocean: Flow patterns, tracers and transports. Progress in Oceanography,Vol. 33, Pergamon, 1–92
- Levitus, S. (1982) Climatological Atlas of the World Ocean, NOAA Professional Paper No. 13, U.S. Govt. Printing Office, 173 pp., -ftp://ftp.nodc.noaa.gov/pub/data.nodc/woa/PUBLICATIONS/levitus_atlas_1982.pdf
- Levitus, S. (1982) Climatological Atlas of the World Ocean, NOAA Professional Paper No. 13, U.S. Govt. Printing Office, 173 pp. - ftp://ftp.nodc.noaa.gov/pub/data.nodc/woa/PUBLICATIONS/levitus_atlas_1982.pdf[permanent dead link]
- Jackett, David R., Trevor J. McDougall, 1997: A Neutral Density Variable for the World's Oceans. J. Phys. Oceanogr., 27, page 248
- Jackett, David R., Trevor J. McDougall, 1997: A Neutral Density Variable for the World's Oceans. J. Phys. Oceanogr., 27, page 252
- Jackett, David R., Trevor J. McDougall, 1997: A Neutral Density Variable for the World's Oceans. J. Phys. Oceanogr., 27, 237–263.
- World Climate Research Programme (WOCW), International Newsletter, June 1995.
- Andreas Klocker, Trevor J. McDougall, David R. Jackett, 2007, “Diapycnal motion due to neutral helicity”).
- Rui Xin Huang, 2010: Is the neutral surface really neutral?
- NOAA, U.S. Department of Commerce, 1982: Climatological Atlas of the World Ocean,ftp://ftp.nodc.noaa.gov/pub/data.nodc/woa/PUBLICATIONS/levitus_atlas_1982.pdf[permanent dead link]