Girinskii potential conversion functions
Revision as of 13:28, 7 December 2007 by Analytophile (talk | contribs) (Added code and description)
- Language
- c++
- Purpose
- Convert between single layer freshwater/saltwater head and potential. From pages 49, 51, and 100 in Strack, 1989
<math>\Phi=\tfrac{1}{2}k\frac{\rho}{\rho_s-\rho}\left[h-\frac{\rho_s}{\rho}H_{sea}+\frac{\rho_s-\rho}{\rho}H\right]^2</math>
<math>h=\frac{\rho_s}{\rho}H_{sea}-\frac{\rho_s-\rho}{\rho}H+\sqrt{\frac{2}{k}\frac{\rho_s-\rho}{\rho}\Phi}</math>
Where <math>\alpha=\frac{\rho_s}{\rho}</math>, and <math>\rho_s</math> and <math>\rho</math> are the density of saltwater and freshwater, respectively. <math>H_{sea}</math> is the sea level elevation and <math>H</math> is the confined layer thickness.
- Author
- James R. Craig
#include <math.h> /******************************************************************************** IsConfined Returns true if local layer condition is confined, false otehrwise -------------------------------------------------------------------------------*/ inline bool IsConfined (const double &potential, const double &K, const double &H){ return (potential>=0.5*K*H*H); }
#include <math.h> /******************************************************************************** ConvertToHead Converts freshwater potential to LOCAL head (head-base) alpha is the ratio of saltwater density to freshwater ~1.035 (1.0 for non-intrusion problems) Hsea is the sea level in local head (<=0.0 for non-intrusion problems) -------------------------------------------------------------------------------*/ inline double ConvertToHead (const double &potential, const double &K, const double &H, const double &Hsea, const double &alpha) { if (K==0){return 0.0;} double Hs=max(Hsea,0.0); if (potential>=0.5*K*H*H){ //confined-------------------------- if (potential>=K*H*(alpha*Hs-0.5*H)){ //freshwater head return (potential/(K*H)+ 0.5*H); } else if (potential<=alpha*K*H*(Hs-0.5*H)){ //saltwater head return alpha*Hs-(alpha-1.0)*H; } else { //mixed head return sqrt(2.0*alpha*(potential-alpha*K*H*(Hs-0.5*H))/K*(1-alpha)); } } else { //unconfined------------------------- if (potential>0.5*K*alpha*alpha*Hs*Hs){ //freshwater head return sqrt(2.0*max(potential,0.0)/K); } else if (potential<0.5*K*alpha*Hs*Hs){ //saltwater head return Hs; } else{ //mixed head return sqrt(2.0*(alpha-1.0)*(max(potential,0.0)-0.5*K*alpha*Hs*Hs)/(K*alpha))+Hs; } } }
#include <math.h> /******************************************************************************** ConvertToPotential Converts from LOCAL freshwater head to potential alpha is the ratio of saltwater density to freshwater ~1.035 (1.0 for non-intrusion problems) Hsea is the sea level in local head (<=0.0 for non-intrusion problems) -------------------------------------------------------------------------------*/ inline double ConvertToPotential(const double &Head, const double &K, const double &H, const double &Hsea, const double &alpha) { double Hs=max(Hsea,0.0); if (Head>H){ //confined-------------------------- if (Head>=alpha*Hs){ //freshwater Head return K*H*(Head-0.5*H); } else{ //saltwater Head return 0.5*K*(1-alpha)/alpha*pow(Head+(alpha-1)*H-alpha*Hs,2)+alpha*K*H*(Hs-0.5*H); } } else{ //unconfined------------------------- if (Head>=alpha*Hs){ //freshwater Head return 0.5*K*Head*Head; } else{ //saltwater Head return 0.5*K*(alpha/(alpha-1.0))*pow((Head-Hs),2)+0.5*K*alpha*Hs*Hs; } } }