Changeset 196 for codes/icosagcm
- Timestamp:
- 07/06/14 09:10:13 (10 years ago)
- Location:
- codes/icosagcm/trunk/src
- Files:
-
- 1 added
- 3 edited
Legend:
- Unmodified
- Added
- Removed
-
codes/icosagcm/trunk/src/physics.f90
r186 r196 1 1 MODULE physics_mod 2 3 USE field_mod 4 5 PRIVATE 6 7 INTEGER, PARAMETER :: phys_none=0, phys_HS94=1, phys_DCMIP=2 8 9 INTEGER :: nb_extra_physics_2D, nb_extra_physics_3D, phys_type 10 TYPE(t_field),POINTER :: f_extra_physics_2D(:), f_extra_physics_3D(:) 2 11 3 12 CHARACTER(LEN=255) :: physics_type="automatic" 4 13 !$OMP THREADPRIVATE(physics_type) 5 14 15 PUBLIC :: physics, init_physics 6 16 CONTAINS 7 17 … … 10 20 USE etat0_mod 11 21 USE icosa 12 USE physics_ dcmip_mod,init_physics_dcmip=>init_physics13 ! USE physics_dry_mod 22 USE physics_interface_mod 23 USE physics_dcmip_mod, ONLY : init_physics_dcmip=>init_physics 14 24 IMPLICIT NONE 15 25 … … 23 33 SELECT CASE(TRIM(etat0_type)) 24 34 CASE('held_suarez') 35 phys_type = phys_HS94 25 36 CASE DEFAULT 26 37 IF(is_mpi_root) PRINT*,"NO PHYSICAL PACKAGE USED" 38 phys_type = phys_none 27 39 END SELECT 28 29 40 30 41 CASE ('dcmip') 31 42 CALL init_physics_dcmip 32 33 CASE ('dry') 34 ! CALL init_physics_dry 35 43 nb_extra_physics_2D=1 ! precl 44 nb_extra_physics_3D=0 45 phys_type = phys_DCMIP 36 46 CASE DEFAULT 37 PRINT*, 'init_physics : Bad selector for variable physics <',TRIM(physics_type),&38 '> options are <automatic>, <dcmip>, <dry>'47 IF(is_mpi_root) PRINT*, 'init_physics : Bad selector for variable physics <',& 48 TRIM(physics_type), '> options are <automatic>, <dcmip>, <dry>' 39 49 STOP 40 50 END SELECT 41 51 52 IF(is_mpi_root) THEN 53 PRINT *, 'phys_type = ',phys_type 54 PRINT *, 'nb_extra_physics_2D = ', nb_extra_physics_2D 55 PRINT *, 'nb_extra_physics_3D = ', nb_extra_physics_3D 56 END IF 57 IF(nb_extra_physics_2D>0) CALL allocate_field(f_extra_physics_2D,field_t,type_real,nb_extra_physics_2D) 58 IF(nb_extra_physics_3D>0) CALL allocate_field(f_extra_physics_3D,field_t,type_real,llm,nb_extra_physics_3D) 59 42 60 END SUBROUTINE init_physics 43 61 44 62 SUBROUTINE physics(it,jD_cur,jH_cur,f_phis, f_ps, f_theta_rhodz, f_ue, f_q) 45 63 USE icosa 46 ! USE physics_dry_mod 47 USE physics_dcmip_mod, physics_dcmip=>physics 48 USE etat0_mod 64 USE physics_interface_mod 65 ! USE etat0_mod 49 66 USE etat0_heldsz_mod 50 67 IMPLICIT NONE … … 56 73 TYPE(t_field),POINTER :: f_ue(:) 57 74 TYPE(t_field),POINTER :: f_q(:) 75 REAL(rstd),POINTER :: phis(:) 76 REAL(rstd),POINTER :: ps(:) 77 REAL(rstd),POINTER :: theta_rhodz(:,:) 78 REAL(rstd),POINTER :: ue(:,:) 79 REAL(rstd),POINTER :: q(:,:,:) 80 58 81 LOGICAL:: firstcall,lastcall 59 60 SELECT CASE(TRIM(physics_type)) 61 CASE ('automatic') 62 63 SELECT CASE(TRIM(etat0_type)) 64 CASE('held_suarez') 65 ! CALL transfert_request(f_ps,req_i1) 66 ! CALL transfert_request(f_theta_rhodz,req_i1) 67 ! CALL transfert_request(f_ue,req_e1_vect) 68 CALL held_suarez(f_ps,f_theta_rhodz,f_ue) 69 CASE DEFAULT 70 ! PRINT*,"NO PHYSICAL PACAKAGE USED" ! FIXME MPI 71 END SELECT 72 73 CASE ('dcmip') 74 CALL physics_dcmip(it, f_phis, f_ps, f_theta_rhodz, f_ue, f_q) 75 76 CASE ('dry') 77 ! CALL physics_dry(it,jD_cur,jH_cur,f_phis, f_ps, f_theta_rhodz, f_ue, f_q) 78 82 INTEGER :: ind 83 TYPE(physics_inout) :: args 84 85 SELECT CASE(phys_type) 86 CASE (phys_none) 87 ! No physics, do nothing 88 CASE(phys_HS94) 89 CALL held_suarez(f_ps,f_theta_rhodz,f_ue) 79 90 CASE DEFAULT 80 PRINT*, 'Bad selector for variable physics <',TRIM(physics_type), & 81 '> options are <automatic>, <dcmip>, <dry>' 82 STOP 91 CALL transfert_request(f_ps,req_i1) 92 CALL transfert_request(f_theta_rhodz,req_i1) 93 CALL transfert_request(f_ue,req_e1_vect) 94 CALL transfert_request(f_q,req_i1) 95 96 args%dt_phys = dt 97 DO ind=1,ndomain 98 IF (.NOT. assigned_domain(ind)) CYCLE 99 CALL swap_dimensions(ind) 100 CALL swap_geometry(ind) 101 102 phis=f_phis(ind) 103 ps=f_ps(ind) 104 theta_rhodz=f_theta_rhodz(ind) 105 ue=f_ue(ind) 106 q=f_q(ind) 107 108 IF(nb_extra_physics_2D>0) args%extra_2D=f_extra_physics_2D(ind) 109 IF(nb_extra_physics_3D>0) args%extra_3D=f_extra_physics_3D(ind) 110 CALL physics_column(args, phis, ps, theta_rhodz, ue, q) 111 ENDDO 112 113 IF (mod(it,itau_out)==0 ) THEN 114 IF(nb_extra_physics_2D>0) CALL writefield("extra_physics_2D",f_extra_physics_2D) 115 IF(nb_extra_physics_3D>0) CALL writefield("extra_physics_3D",f_extra_physics_3D) 116 ENDIF 83 117 END SELECT 84 118 85 119 END SUBROUTINE physics 86 120 121 SUBROUTINE physics_column(args, phis, ps, theta_rhodz, ue, q) 122 USE icosa 123 USE wind_mod 124 USE pression_mod 125 USE theta2theta_rhodz_mod 126 USE physics_interface_mod 127 USE physics_dcmip_mod 128 IMPLICIT NONE 129 TYPE(physics_inout) :: args 130 REAL(rstd) :: phis(iim*jjm) 131 REAL(rstd) :: ps(iim*jjm) 132 REAL(rstd) :: theta_rhodz(iim*jjm,llm) 133 REAL(rstd) :: ue(3*iim*jjm,llm) 134 REAL(rstd), TARGET :: q(iim*jjm,llm,nqtot) 135 ! local arrays 136 REAL(rstd), TARGET :: lat(iim*jjm) 137 REAL(rstd), TARGET :: lon(iim*jjm) 138 REAL(rstd), TARGET :: p(iim*jjm,llm+1) 139 REAL(rstd), TARGET :: Temp(iim*jjm,llm) 140 REAL(rstd), TARGET :: ulon(iim*jjm,llm) 141 REAL(rstd), TARGET :: ulat(iim*jjm,llm) 142 REAL(rstd), TARGET :: dTemp(iim*jjm,llm) 143 REAL(rstd), TARGET :: dulon(iim*jjm,llm) 144 REAL(rstd), TARGET :: dulat(iim*jjm,llm) 145 REAL(rstd), TARGET :: dq(iim*jjm,llm,nqtot) 146 REAL(rstd) :: uc(iim*jjm,3,llm) ! 3D velocity at cell centers 147 148 INTEGER :: i,j,ij,l 149 REAL(rstd) :: due, dt2 150 151 DO j=jj_begin,jj_end 152 DO i=ii_begin,ii_end 153 ij=(j-1)*iim+i 154 CALL xyz2lonlat(xyz_i(ij,:),lon(ij),lat(ij)) 155 ENDDO 156 ENDDO 157 158 ! Reconstruct wind vector at hexagons 159 CALL compute_pression(ps,p,0) 160 CALL compute_theta_rhodz2temperature(ps,theta_rhodz,Temp,0) 161 CALL compute_wind_centered(ue,uc) 162 CALL compute_wind_centered_lonlat_compound(uc, ulon, ulat) 163 args%ngrid = iim*jjm 164 args%lon => lon 165 args%lat => lat 166 args%p => p 167 args%Temp => Temp 168 args%ulon => ulon 169 args%ulat => ulat 170 args%q => q 171 args%dTemp => dTemp 172 args%dulon => dulon 173 args%dulat => dulat 174 args%dq => dq 175 176 SELECT CASE(phys_type) 177 CASE (phys_DCMIP) 178 CALL compute_phys_wrap(args) 179 END SELECT 180 181 q = q + args%dt_phys * dq 182 Temp = Temp + args%dt_phys * dTemp 183 CALL compute_temperature2theta_rhodz(ps,Temp,theta_rhodz,0) 184 185 ! Reconstruct wind tendencies at edges and add 186 CALL compute_wind_centered_from_lonlat_compound(dulon,dulat,uc) 187 dt2=.5*args%dt_phys 188 DO l=1,llm 189 DO j=jj_begin,jj_end 190 DO i=ii_begin,ii_end 191 ij=(j-1)*iim+i 192 due = sum((uc(ij,:,l) + uc(ij+t_right,:,l))*ep_e(ij+u_right,:) ) 193 ue(ij+u_right,l) = ue(ij+u_right,l) + dt2*due 194 195 due = sum( (uc(ij,:,l) + uc(ij+t_lup,:,l))*ep_e(ij+u_lup,:) ) 196 ue(ij+u_lup,l)=ue(ij+u_lup,l) + dt2*due 197 198 due = sum( (uc(ij,:,l) + uc(ij+t_ldown,:,l))*ep_e(ij+u_ldown,:) ) 199 ue(ij+u_ldown,l)=ue(ij+u_ldown,l) + dt2*due 200 ENDDO 201 ENDDO 202 ENDDO 203 204 END SUBROUTINE physics_column 205 87 206 END MODULE physics_mod -
codes/icosagcm/trunk/src/physics_dcmip.f90
r192 r196 6 6 !$OMP THREADPRIVATE(testcase) 7 7 8 PUBLIC init_physics, physics9 8 TYPE(t_field),POINTER :: f_out_i(:) 10 9 TYPE(t_field),POINTER :: f_precl(:) 11 REAL(rstd),POINTER :: out_i(:,:) 10 REAL(rstd),POINTER :: out_i(:,:) 11 12 PUBLIC :: compute_phys_wrap, init_physics 12 13 13 14 CONTAINS 14 15 15 16 SUBROUTINE init_physics 16 IMPLICIT NONE 17 17 IMPLICIT NONE 18 18 testcase=1 ! OK for 4.2 (moist baroclinic instability) 19 19 CALL getin("dcmip_physics",testcase) 20 CALL allocate_field(f_out_i,field_t,type_real,llm)21 CALL allocate_field(f_precl,field_t,type_real)22 23 20 END SUBROUTINE init_physics 24 21 25 26 SUBROUTINE physics( it, f_phis, f_ps, f_theta_rhodz, f_ue, f_q) 27 USE icosa 28 IMPLICIT NONE 29 INTEGER,INTENT(IN) :: it 30 TYPE(t_field),POINTER :: f_phis(:) 31 TYPE(t_field),POINTER :: f_ps(:) 32 TYPE(t_field),POINTER :: f_theta_rhodz(:) 33 TYPE(t_field),POINTER :: f_ue(:) 34 TYPE(t_field),POINTER :: f_q(:) 35 36 REAL(rstd),POINTER :: phis(:) 37 REAL(rstd),POINTER :: ps(:) 38 REAL(rstd),POINTER :: theta_rhodz(:,:) 39 REAL(rstd),POINTER :: ue(:,:) 40 REAL(rstd),POINTER :: q(:,:,:) 41 REAL(rstd),POINTER :: precl(:) 42 INTEGER :: ind 43 44 CALL transfert_request(f_ue,req_e1_vect) 45 DO ind=1,ndomain 46 IF (.NOT. assigned_domain(ind)) CYCLE 47 CALL swap_dimensions(ind) 48 CALL swap_geometry(ind) 49 50 phis=f_phis(ind) 51 ps=f_ps(ind) 52 theta_rhodz=f_theta_rhodz(ind) 53 ue=f_ue(ind) 54 q=f_q(ind) 55 out_i=f_out_i(ind) 56 precl=f_precl(ind) 57 58 CALL compute_physics( phis, ps, theta_rhodz, ue, q(:,:,1), precl) 59 60 ENDDO 61 62 ! CALL writefield("out_i",f_out_i) 63 64 IF (mod(it,itau_out)==0 ) THEN 65 CALL writefield("precl",f_precl) 66 ENDIF 67 68 END SUBROUTINE physics 69 70 SUBROUTINE compute_physics(phis, ps, theta_rhodz, ue, q, precl) 71 USE icosa 72 USE pression_mod 73 USE exner_mod 74 USE theta2theta_rhodz_mod 75 USE geopotential_mod 76 USE wind_mod 77 IMPLICIT NONE 78 REAL(rstd) :: phis(iim*jjm) 79 REAL(rstd) :: ps(iim*jjm) 80 REAL(rstd) :: theta_rhodz(iim*jjm,llm) 81 REAL(rstd) :: ue(3*iim*jjm,llm) 82 REAL(rstd) :: q(iim*jjm,llm) 83 REAL(rstd) :: precl(iim*jjm) 84 85 REAL(rstd) :: p(iim*jjm,llm+1) 86 REAL(rstd) :: pks(iim*jjm) 87 REAL(rstd) :: pk(iim*jjm,llm) 88 REAL(rstd) :: phi(iim*jjm,llm) 89 REAL(rstd) :: T(iim*jjm,llm) 90 REAL(rstd) :: Tfi(iim*jjm,llm) 91 REAL(rstd) :: theta(iim*jjm,llm) 92 93 REAL(rstd) :: uc(iim*jjm,3,llm) 94 REAL(rstd) :: u(iim*jjm,llm) 95 REAL(rstd) :: v(iim*jjm,llm) 96 REAL(rstd) :: ufi(iim*jjm,llm) 97 REAL(rstd) :: vfi(iim*jjm,llm) 98 REAL(rstd) :: qfi(iim*jjm,llm) 99 REAL(rstd) :: utemp(iim*jjm,llm) 100 REAL(rstd) :: vtemp(iim*jjm,llm) 101 REAL(rstd) :: lat(iim*jjm) 102 REAL(rstd) :: lon 103 REAL(rstd) :: pmid(iim*jjm,llm) 104 REAL(rstd) :: pint(iim*jjm,llm+1) 105 REAL(rstd) :: pdel(iim*jjm,llm) 106 INTEGER :: i,j,l,ij 107 108 ! PRINT *,'Entering in DCMIP physics' 109 CALL compute_pression(ps,p,0) 110 CALL compute_exner(ps,p,pks,pk,0) 111 CALL compute_theta_rhodz2theta(ps,theta_rhodz,theta,0) 112 CALL compute_geopotential(phis,pks,pk,theta,phi,0) 113 CALL compute_theta_rhodz2temperature(ps,theta_rhodz,T,0) 114 ! Reconstruct wind vector at hexagons 115 CALL compute_wind_centered(ue,uc) 116 CALL compute_wind_centered_lonlat_compound(uc, u, v) 117 118 ! Convert from Tv to T 119 DO l=1,llm 120 DO j=jj_begin,jj_end 121 DO i=ii_begin,ii_end 122 ij=(j-1)*iim+i 123 T(ij,l)=T(ij,l)/(1+0.608*q(ij,l)) 124 ENDDO 125 ENDDO 126 ENDDO 127 128 DO j=jj_begin,jj_end 129 DO i=ii_begin,ii_end 130 ij=(j-1)*iim+i 131 CALL xyz2lonlat(xyz_i(ij,:),lon,lat(ij)) 132 ENDDO 133 ENDDO 134 135 ! bottom-up indexing (DYNAMICO) : u,utemp, v,vtemp 136 ! top-down vertical indexing (DCMIP) : ufi, vfi, ... 137 ! => copy fields and mirror vertical index 22 SUBROUTINE compute_phys_wrap(args) 23 USE physics_interface_mod 24 TYPE(physics_inout) :: args 25 CALL compute_physics(args%ngrid, args%dt_phys, args%lat, & 26 args%p, args%Temp, args%ulon, args%ulat, args%q(:,:,1), & 27 args%dTemp, args%dulon, args%dulat, args%dq(:,:,1), args%extra_2D(:,1)) 28 END SUBROUTINE compute_phys_wrap 29 30 SUBROUTINE compute_physics(ngrid,dt_phys,lat, p,Temp,u,v,q, dTemp,du,dv,dq, precl) 31 USE icosa 32 IMPLICIT NONE 33 INTEGER :: ngrid 34 REAL(rstd) :: lat(ngrid) 35 REAL(rstd) :: ps(ngrid) 36 REAL(rstd) :: precl(ngrid) 37 ! arguments with bottom-up indexing (DYNAMICO) 38 REAL(rstd) :: p(ngrid,llm+1) 39 REAL(rstd) :: Temp(ngrid,llm) 40 REAL(rstd) :: u(ngrid,llm) 41 REAL(rstd) :: v(ngrid,llm) 42 REAL(rstd) :: q(ngrid,llm) 43 REAL(rstd) :: dTemp(ngrid,llm) 44 REAL(rstd) :: du(ngrid,llm) 45 REAL(rstd) :: dv(ngrid,llm) 46 REAL(rstd) :: dq(ngrid,llm) 47 ! local arrays with top-down vertical indexing (DCMIP) 48 REAL(rstd) :: pint(ngrid,llm+1) 49 REAL(rstd) :: pmid(ngrid,llm) 50 REAL(rstd) :: pdel(ngrid,llm) 51 REAL(rstd) :: Tfi(ngrid,llm) 52 REAL(rstd) :: ufi(ngrid,llm) 53 REAL(rstd) :: vfi(ngrid,llm) 54 REAL(rstd) :: qfi(ngrid,llm) 55 56 INTEGER :: i,j,l,ll,ij 57 REAL(rstd) :: dt_phys, inv_dt 58 59 ! prepare input fields and mirror vertical index 60 ps(:) = p(:,1) ! surface pressure 61 138 62 DO l=1,llm+1 139 63 DO j=jj_begin,jj_end … … 144 68 ENDDO 145 69 ENDDO 146 147 ! Pressure inside layers 70 148 71 DO l=1,llm 149 DO j=jj_begin,jj_end 150 DO i=ii_begin,ii_end 151 ij=(j-1)*iim+i 152 pmid(ij,l)=0.5*(pint(ij,l)+pint(ij,l+1)) 153 ENDDO 154 ENDDO 72 ll=llm+1-l 73 DO j=jj_begin,jj_end 74 DO i=ii_begin,ii_end 75 ij=(j-1)*iim+i 76 pmid(ij,l)=0.5*(pint(ij,l)+pint(ij,l+1)) ! Pressure inside layers 77 pdel(ij,l)=pint(ij,l+1)-pint(ij,l) ! Pressure difference between two layers 78 ufi(ij,l)=u(ij,ll) 79 vfi(ij,l)=v(ij,ll) 80 qfi(ij,l)=q(ij,ll) 81 Tfi(ij,l)=Temp(ij,ll)/(1+0.608*qfi(ij,l)) 82 ENDDO 83 ENDDO 155 84 ENDDO 156 85 157 ! Pressure difference between two layers 86 CALL simple_physics(ngrid, llm, dt_phys, lat, tfi, qfi , ufi, vfi, pmid, pint, pdel, 1./pdel, ps, precl, testcase) 87 88 ! Mirror vertical index and compute tendencies 89 inv_dt = 1./dt_phys 158 90 DO l=1,llm 159 DO j=jj_begin,jj_end 160 DO i=ii_begin,ii_end 161 ij=(j-1)*iim+i 162 pdel(ij,l)=pint(ij,l+1)-pint(ij,l) 163 ENDDO 164 ENDDO 165 ENDDO 166 167 ! Copy T,u,v,q for input to physics 168 DO l=1,llm 169 DO j=jj_begin,jj_end 170 DO i=ii_begin,ii_end 171 ij=(j-1)*iim+i 172 Tfi(ij,l)=T(ij,llm+1-l) 173 ufi(ij,l)=u(ij,llm+1-l) 174 vfi(ij,l)=v(ij,llm+1-l) 175 qfi(ij,l)=q(ij,llm+1-l) 176 ENDDO 177 ENDDO 178 ENDDO 179 180 CALL simple_physics(iim*jjm, llm, dt, lat, tfi, qfi , ufi, vfi, pmid, pint, pdel, 1/pdel, ps, precl, testcase) 181 182 ! Copy back T,u,v,q and mirror vertical index 183 DO l=1,llm 184 DO j=jj_begin,jj_end 185 DO i=ii_begin,ii_end 186 ij=(j-1)*iim+i 187 T(ij,l)=Tfi(ij,llm+1-l) 188 utemp(ij,l)=ufi(ij,llm+1-l) 189 vtemp(ij,l)=vfi(ij,llm+1-l) 190 q(ij,l)=qfi(ij,llm+1-l) 191 ENDDO 192 ENDDO 193 ENDDO 194 195 196 ! Convert back T to Tv 197 DO l=1,llm 198 DO j=jj_begin,jj_end 199 DO i=ii_begin,ii_end 200 ij=(j-1)*iim+i 201 T(ij,l)=T(ij,l)*(1+0.608*q(ij,l)) 202 ENDDO 203 ENDDO 204 ENDDO 205 206 ! Compute velocity update at hexagons 207 utemp=utemp-u 208 vtemp=vtemp-v 209 210 ! lon-lat -> 3D 211 DO l=1,llm 212 DO j=jj_begin,jj_end 213 DO i=ii_begin,ii_end 214 ij=(j-1)*iim+i 215 uc(ij,:,l)=utemp(ij,l)*elon_i(ij,:)+vtemp(ij,l)*elat_i(ij,:) 216 ENDDO 217 ENDDO 91 ll=llm+1-l 92 DO j=jj_begin,jj_end 93 DO i=ii_begin,ii_end 94 ij=(j-1)*iim+i 95 dTemp(ij,l) = inv_dt * ( Tfi(ij,ll)*(1+0.608*qfi(ij,ll)) - Temp(ij,l) ) 96 du(ij,l) = inv_dt * (ufi(ij,ll) - u(ij,l)) 97 dv(ij,l) = inv_dt * (vfi(ij,ll) - v(ij,l)) 98 dq(ij,l) = inv_dt * (qfi(ij,ll) - q(ij,l)) 99 ENDDO 100 ENDDO 218 101 ENDDO 219 220 ! Update velocity at velocity points221 DO l=1,llm222 DO j=jj_begin,jj_end223 DO i=ii_begin,ii_end224 ij=(j-1)*iim+i225 ue(ij+u_right,l)=ue(ij+u_right,l)+sum( 0.5*(uc(ij,:,l) + uc(ij+t_right,:,l))*ep_e(ij+u_right,:) )226 ue(ij+u_lup,l)=ue(ij+u_lup,l)+sum( 0.5*(uc(ij,:,l) + uc(ij+t_lup,:,l))*ep_e(ij+u_lup,:) )227 ue(ij+u_ldown,l)=ue(ij+u_ldown,l)+sum( 0.5*(uc(ij,:,l) + uc(ij+t_ldown,:,l))*ep_e(ij+u_ldown,:) )228 ENDDO229 ENDDO230 ENDDO231 232 CALL compute_temperature2theta_rhodz(ps,T,theta_rhodz,0)233 234 102 235 103 END SUBROUTINE compute_physics 236 104 105 END MODULE physics_dcmip_mod 237 106 238 107 subroutine simple_physics (pcols, pver, dtime, lat, t, q, u, v, pmid, pint, pdel, rpdel, ps, precl, test) … … 501 370 do i=1,pcols 502 371 qsat = epsilo*e0/pmid(i,k)*exp(-latvap/rh2o*((1._r8/t(i,k))-1._r8/T0)) ! saturation specific humidity 503 out_i(i,llm+1-k)=q(i,k)-qsat372 ! out_i(i,llm+1-k)=q(i,k)-qsat 504 373 if (q(i,k) > qsat) then ! saturated? 505 374 tmp = 1._r8/dtime*(q(i,k)-qsat)/(1._r8+(latvap/cpair)*(epsilo*latvap*qsat/(rair*t(i,k)**2))) … … 694 563 end subroutine simple_physics 695 564 696 697 698 699 END MODULE physics_dcmip_mod -
codes/icosagcm/trunk/src/wind.f90
r186 r196 185 185 END SUBROUTINE compute_wind_from_lonlat_compound 186 186 187 SUBROUTINE compute_wind_centered_from_lonlat_compound(ulon, ulat, u) 188 USE icosa 189 190 IMPLICIT NONE 191 REAL(rstd) :: u(iim*jjm,3,llm) 192 REAL(rstd) :: ulon(iim*jjm,llm) 193 REAL(rstd) :: ulat(iim*jjm,llm) 194 195 INTEGER :: i,j,ij,l 196 DO l=1,llm 197 DO j=jj_begin-1,jj_end+1 198 DO i=ii_begin-1,ii_end+1 199 ij=(j-1)*iim+i 200 u(ij,:,l)=ulon(ij,l)*elon_i(ij,:)+ ulat(ij,l)*elat_i(ij,:) 201 ENDDO 202 ENDDO 203 ENDDO 204 205 END SUBROUTINE compute_wind_centered_from_lonlat_compound 206 187 207 SUBROUTINE compute_wind2D_from_lonlat_compound(ulon, ulat, u) 188 208 USE icosa
Note: See TracChangeset
for help on using the changeset viewer.