Changeset 67 for codes/icosagcm/trunk
- Timestamp:
- 08/01/12 06:58:27 (12 years ago)
- Location:
- codes/icosagcm/trunk/src
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
codes/icosagcm/trunk/src/etat0.f90
r62 r67 37 37 CASE DEFAULT 38 38 PRINT*, 'Bad selector for variable etat0 <',etat0_type, & 39 '> options are <jablonowsky06>, <academic>, <dcmip[1- 3]> '39 '> options are <jablonowsky06>, <academic>, <dcmip[1-4]> ' 40 40 STOP 41 41 END SELECT -
codes/icosagcm/trunk/src/etat0_dcmip41.f90
r62 r67 82 82 83 83 DO l=1,llm 84 eta(l)= 0.5 *( ap(l)/p s0+bp(l) + ap(l+1)/ps0+bp(l+1) )84 eta(l)= 0.5 *( ap(l)/preff+bp(l) + ap(l+1)/preff+bp(l+1) ) 85 85 etav(l)=(eta(l)-eta0)*Pi/2 86 86 ENDDO 87 etas=ap(1) +bp(1)87 etas=ap(1)/preff+bp(1) 88 88 etavs=(etas-eta0)*Pi/2 89 89 … … 154 154 phis(ij)=phis_ave+u0*cos(etavs)**1.5*( (-2*sin(lat)**6 * (cos(lat)**2+1./3) + 10./63 )*u0*cos(etavs)**1.5 & 155 155 +(8./5*cos(lat)**3 * (sin(lat)**2 + 2./3) - Pi/4)*radius*Omega ) 156 ! phis(ij)=phis_ave+u0*cos(etavs)**1.5 157 156 158 ENDDO 157 159 ENDDO … … 170 172 ij=(j-1)*iim+i 171 173 CALL xyz2lonlat(xyz_i(ij,:)/radius,lon,lat) 172 dthetaodeta=dthetaodeta_ave + 0.25 * Pi * u0/Rd*(1-kappa)*eta(l)**(-kappa)*sin(etav(l))*cos(etav(l))**0.5 * Y(ij,l) & 173 - 3./16. * Pi**2 * u0 /Rd * eta(l)**(1-kappa) * sin(etav(l))**2 * cos(etav(l))**0.5 *Y(ij,l)& 174 - 9./8. * Pi**2 * u0 /Rd * eta(l)**(1-kappa) * sin(etav(l))**2 * cos(etav(l)) & 174 dthetaodeta=dthetaodeta_ave + 3./4. * Pi * u0/Rd*(1-kappa)*eta(l)**(-kappa)*sin(etav(l))*cos(etav(l))**0.5 * Y(ij,l) & 175 + 3/8. * Pi**2*u0/Rd * eta(l)**(1-kappa) * cos(etav(l))**1.5 * Y(ij,l) & 176 - 3./16. * Pi**2 * u0 /Rd * eta(l)**(1-kappa) * sin(etav(l))**2 * cos(etav(l))**(-0.5) *Y(ij,l)& 177 - 9./8. * Pi**2 * u0 /Rd * eta(l)**(1-kappa) * sin(etav(l))**2 * cos(etav(l)) & 175 178 * (-2*sin(lat)**6*(cos(lat)**2+1./3.)+10./63.) 176 179 dthetaodlat=3./4.*Pi*u0/Rd*eta(l)**(1-kappa)*sin(etav(l))*cos(etav(l))**0.5 & … … 182 185 K=sin(latc)*sin(lat)+cos(latc)*cos(lat)*cos(lon-lonc) 183 186 r=radius*acos(K) 184 vort = -4*u0/radius*cos(etav(l))**1.5 * sin(lat) * cos(lat) * (2 -5*sin(lat)**2) &185 + up0/radius*exp(-(r/(0.1*radius))**2) * (tan(lat)-2*(r /(0.1*radius))**2 * acos(K) * (sin(latc)*cos(lat) &187 vort = -4*u0/radius*cos(etav(l))**1.5 * sin(lat) * cos(lat) * (2.-5.*sin(lat)**2) & 188 + up0/radius*exp(-(r/(0.1*radius))**2) * (tan(lat)-2*(radius/(0.1*radius))**2 * acos(K) * (sin(latc)*cos(lat) & 186 189 -cos(latc)*sin(lat)*cos(lon-lonc))/(sqrt(1-K**2))) 187 duodeta= -u0*sin(2*lat)*3*Pi/4*cos(etav(l))**0.5*sin(etav(l)) 188 189 q(ij,l,2)=g/preff*(-1./radius*duodeta*dthetaodlat-(2*sin(lat)*omega+vort)*dthetaodeta) 190 q(ij,l,2)=ABS(g/preff*(-1./radius*duodeta*dthetaodlat-(2*sin(lat)*omega+vort)*dthetaodeta)) 190 191 ENDDO 191 192 ENDDO
Note: See TracChangeset
for help on using the changeset viewer.