| 1 |
function [qsat,dQdT]=calc_Qsat(kpot,p,t); |
| 2 |
% [qsat,dQdT]=calc_Qsat(kpot,p,t); |
| 3 |
% kpot=0 : t = T = absolute Temp ; kpot=1 : t = theta = pot.temp |
| 4 |
% p = normalized pressure = p/1000.mb |
| 5 |
% output: qsat (g/kg) ; d.qsat/d.T |
| 6 |
% ktyp=1 : p in mb ; t = absolute Temp ; |
| 7 |
|
| 8 |
g=9.81; |
| 9 |
kappa=2./7.; |
| 10 |
|
| 11 |
np=prod(size(p)); |
| 12 |
nt=prod(size(t)); |
| 13 |
nx=nt/np; |
| 14 |
if nx ~= fix(nx) | nx < 1 |
| 15 |
fprintf(' Pb in dimensions: dim_p,dim_t= %i , %i ; ratio= %f \n',np,nt,nx); |
| 16 |
qsat=0.; dQdT=0.; |
| 17 |
return |
| 18 |
end |
| 19 |
|
| 20 |
T=reshape(t,nx,np); |
| 21 |
P=reshape(p,1,np); |
| 22 |
if kpot == 1, |
| 23 |
factp=P.^kappa; |
| 24 |
fp=ones(nx,1)*factp; |
| 25 |
size(fp); |
| 26 |
T=fp.*T; |
| 27 |
end |
| 28 |
|
| 29 |
e0=6.108e-3 ; c1=17.269 ; c2=21.875 ; t0c=273.16 ; ct1=35.86 ; ct2=7.66 ; |
| 30 |
c0q=622. ; c1q=0.378 ; |
| 31 |
|
| 32 |
qsat=zeros(nx,np); |
| 33 |
dQdT=zeros(nx,np); |
| 34 |
|
| 35 |
for k=1:np, |
| 36 |
for j=1:nx, |
| 37 |
if T(j,k) >= t0c |
| 38 |
qE=e0*exp(c1*(T(j,k)-t0c)/(T(j,k)-ct1)) ; |
| 39 |
qsat(j,k)=c0q*qE/(P(k)-c1q*qE); |
| 40 |
dQdT(j,k)=c1*P(k)*qsat(j,k)/(P(k)-c1q*qE)*(t0c-ct1)/(T(j,k)-ct1)^2 ; |
| 41 |
elseif T(j,k) > ct2 |
| 42 |
qE=e0*exp(c2*(T(j,k)-t0c)/(T(j,k)-ct2)) ; |
| 43 |
qsat(j,k)=c0q*qE/(P(k)-c1q*qE); |
| 44 |
dQdT(j,k)=c2*P(k)*qsat(j,k)/(P(k)-c1q*qE)*(t0c-ct2)/(T(j,k)-ct2)^2 ; |
| 45 |
end |
| 46 |
end |
| 47 |
end |
| 48 |
|
| 49 |
return |
Parent Directory
| 
Revision Log
| 
Revision Graph