Home > atmlab > h2o > parametrisations > sphere_part_fallspeed.m

sphere_part_fallspeed

PURPOSE ^

SPHERE_PART_FALLSPEED Calculate spherical particle fallspeed in air

SYNOPSIS ^

function [u]=sphere_part_fallspeed(r,T,P)

DESCRIPTION ^

SPHERE_PART_FALLSPEED   Calculate spherical particle fallspeed in air

          Calculates spherical particles fallspeed
          according to R.R. Rogers, "A short course
          in cloud physics"
    
 FORMAT   [u]=sphere_part_fallspeed(r,T,P)

 IN    r   particle radius  (vector) [m]
       T   temperature      (scalar) [K], allowed range is [233.15,303.15]
       P   pressure         (scalar) [Pa]

 OUT   u   fallspeed       [m/s]
         
 example usage: 
               r=[30 50 100 1000]*1e-6;
               T=290;P=1e5;
               [u]=sphere_part_fallspeed(r,T,P)

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

DOWNLOAD ^

sphere_part_fallspeed.m

SOURCE CODE ^

0001 %SPHERE_PART_FALLSPEED   Calculate spherical particle fallspeed in air
0002 %
0003 %          Calculates spherical particles fallspeed
0004 %          according to R.R. Rogers, "A short course
0005 %          in cloud physics"
0006 %
0007 % FORMAT   [u]=sphere_part_fallspeed(r,T,P)
0008 %
0009 % IN    r   particle radius  (vector) [m]
0010 %       T   temperature      (scalar) [K], allowed range is [233.15,303.15]
0011 %       P   pressure         (scalar) [Pa]
0012 %
0013 % OUT   u   fallspeed       [m/s]
0014 %
0015 % example usage:
0016 %               r=[30 50 100 1000]*1e-6;
0017 %               T=290;P=1e5;
0018 %               [u]=sphere_part_fallspeed(r,T,P)
0019 
0020 % 2009-03-30 Created by Bengt Rydberg
0021 
0022 function [u]=sphere_part_fallspeed(r,T,P)
0023 
0024 %some checks
0025 if T<-40+273.15 | T>273.15+30 & min(r)<=40e-6 
0026   error('Temperature out of range!!!T must be inside (233.15-303.15)') 
0027 end
0028 
0029 if length(T)>1
0030   error('T must be a scalar!!!')
0031 end
0032 
0033 if length(P)>1
0034   error('P must be a scalar!!!')
0035 end
0036 
0037 %dynamic viscosity of air (variation with temperature)
0038 muvec  = [1.512 1.564 1.616 1.667 1.717 1.766 1.815 1.862]*1e-5;
0039 Tvec   = [-40:10:30]+273.15;
0040 mu     = interp1(Tvec,muvec,T);
0041 
0042 rhol   = 1e3;    %density of water
0043 Rp     = 287;    %individual gas constant of dry air
0044 rho_air= P/Rp/T; %density of air
0045 rho0=1.2;        %reference density of air
0046 
0047 %loop over particle sizes
0048 for i=1:length(r)
0049 
0050  if r<=40e-6
0051   u(i)=2/9*9.81*rhol/mu*r(i)^2;
0052  elseif r>40e-6 & r<=0.6e-3
0053   u(i)=8e3*r(i);
0054  else
0055   u(i)=2.2e3*sqrt(rho0/rho_air*r(i)*1e2)/1e2;
0056  end
0057 
0058 end

Generated on Mon 15-Sep-2014 13:31:28 by m2html © 2005