- Timestamp:
- 10/20/14 23:42:26 (10 years ago)
- Location:
- codes/icosagcm/trunk/src
- Files:
-
- 12 edited
Legend:
- Unmodified
- Added
- Removed
-
codes/icosagcm/trunk/src/caldyn_gcm.f90
r266 r286 328 328 INTEGER :: ij 329 329 DO ij=ij_begin_ext,ij_end_ext 330 CALL xyz2lonlat(xyz_e(ij+u_right,:),lon,lat) 331 ulon(ij+u_right)=a*omega*cos(lat) 330 ulon(ij+u_right)=a*omega*cos(lat_e(ij+u_right)) 332 331 ulat(ij+u_right)=0 333 332 334 CALL xyz2lonlat(xyz_e(ij+u_lup,:),lon,lat) 335 ulon(ij+u_lup)=a*omega*cos(lat) 333 ulon(ij+u_lup)=a*omega*cos(lat_e(ij+u_lup)) 336 334 ulat(ij+u_lup)=0 337 335 338 CALL xyz2lonlat(xyz_e(ij+u_ldown,:),lon,lat) 339 ulon(ij+u_ldown)=a*omega*cos(lat) 336 ulon(ij+u_ldown)=a*omega*cos(lat_e(ij+u_ldown)) 340 337 ulat(ij+u_ldown)=0 341 338 END DO -
codes/icosagcm/trunk/src/check_conserve.f90
r266 r286 256 256 ij=(j-1)*iim+i 257 257 IF (domain(ind)%own(i,j)) THEN 258 ! CALL xyz2lonlat(xyz_i(ij,:),lon,lat) 259 ! very bad inlining with intel compiler (>14 ?) 260 xyz(:)=xyz_i(ij,:) 261 xyz(:)=xyz(:)/sqrt(sum(xyz(:)**2)) 262 lat=asin(xyz(3)) 263 lon=atan2(xyz(2),xyz(1)) 264 258 lat=lat_i(ij) 265 259 ang=ang+rad*cos(lat)*masse(ij,l)*(ulon(ij,l)+ radomeg*cos(lat)) 266 260 END IF -
codes/icosagcm/trunk/src/dissip_gcm.f90
r250 r286 540 540 541 541 nn = ij+shift_u 542 CALL xyz2lonlat(xyz_e(nn,:),lon,lat)543 542 zcoords = 1 ; hybrid_eta = .FALSE. ! use z instead of p or hyam/hybm 544 CALL test2_schaer_mountain(lon ,lat,p,z,zcoords,hybrid_eta, &543 CALL test2_schaer_mountain(lon_e(nn),lat_e(nn),p,z,zcoords,hybrid_eta, & 545 544 hyam,hybm,shear,ulon,ulat,w,t,phis,ps,rho,q) 546 545 ! u3d = ulon*elon_e(nn,:) + ulat*elat_e(nn,:) -
codes/icosagcm/trunk/src/etat0.f90
r266 r286 158 158 REAL(rstd),INTENT(OUT) :: q(iim*jjm,llm,nqtot) 159 159 160 REAL(rstd) :: lon_i(iim*jjm)161 REAL(rstd) :: lat_i(iim*jjm)162 160 REAL(rstd) :: temp_i(iim*jjm,llm) 163 161 REAL(rstd) :: ulon_i(iim*jjm,llm) 164 162 REAL(rstd) :: ulat_i(iim*jjm,llm) 165 163 166 REAL(rstd) :: lon_e(3*iim*jjm)167 REAL(rstd) :: lat_e(3*iim*jjm)168 164 REAL(rstd) :: ps_e(3*iim*jjm) 169 165 REAL(rstd) :: mass_e(3*iim*jjm,llm) … … 175 171 176 172 INTEGER :: l,i,j,ij 177 178 DO l=1,llm179 DO j=jj_begin-1,jj_end+1180 DO i=ii_begin-1,ii_end+1181 ij=(j-1)*iim+i182 CALL xyz2lonlat(xyz_i(ij,:)/radius,lon_i(ij),lat_i(ij))183 CALL xyz2lonlat(xyz_e(ij+u_right,:)/radius,lon_e(ij+u_right),lat_e(ij+u_right))184 CALL xyz2lonlat(xyz_e(ij+u_lup,:)/radius,lon_e(ij+u_lup),lat_e(ij+u_lup))185 CALL xyz2lonlat(xyz_e(ij+u_ldown,:)/radius,lon_e(ij+u_ldown),lat_e(ij+u_ldown))186 END DO187 END DO188 END DO189 173 190 174 SELECT CASE (TRIM(etat0_type)) -
codes/icosagcm/trunk/src/etat0_academic.f90
r186 r286 109 109 DO i=ii_begin-1,ii_end+1 110 110 ij=(j-1)*iim+i 111 CALL xyz2lonlat(xyz_i(ij,:),lon,lat) 112 ddsin=sin(lat)-sin(pi/20.) 111 ddsin=sin(lat_i(ij))-sin(pi/20.) 113 112 theta(ij,l)=300.*(1+1./18.*(1.-3.*ddsin*ddsin)+thetarappel) 114 113 ENDDO … … 144 143 ij=(j-1)*iim+i 145 144 146 CALL xyz2lonlat(xyz_e(ij+u_right,:),lon,lat)145 lat=lat_e(ij+u_right) 147 146 IF (abs(sin(lat))<1.e-4) lat=1e-4 148 147 x=cos(lat) ; x=x*x ; x=x*x ; x=x*x 149 148 fact(ij+u_right)=(1.-x)/(2.*omega*sin(lat)) 150 149 151 CALL xyz2lonlat(xyz_e(ij+u_lup,:),lon,lat)150 lat=lat_e(ij+u_lup) 152 151 IF (abs(sin(lat))<1.e-4) lat=1e-4 153 152 x=cos(lat) ; x=x*x ; x=x*x ; x=x*x 154 153 fact(ij+u_lup)=(1.-x)/(2.*omega*sin(lat)) 155 154 156 CALL xyz2lonlat(xyz_e(ij+u_ldown,:),lon,lat)155 lat=lat_e(ij+u_ldown) 157 156 IF (abs(sin(lat))<1.e-4) lat=1e-4 158 157 x=cos(lat) ; x=x*x ; x=x*x ; x=x*x 159 158 fact(ij+u_ldown)=(1.-x)/(2.*omega*sin(lat)) 160 161 ! fact(ij+u_right)=-1162 ! fact(ij+u_lup)=-1163 ! fact(ij+u_ldown)=-1164 165 ! CALL xyz2lonlat(xyz_e(ij+u_ldown,:),lon,lat)166 ! ps(ij)=lat167 159 ENDDO 168 160 ENDDO -
codes/icosagcm/trunk/src/etat0_dcmip1.f90
r271 r286 188 188 DO i=ii_begin-1,ii_end+1 189 189 n=(j-1)*iim+i 190 CALL xyz2lonlat(xyz_i(n,:),lon,lat) 191 CALL dist_lonlat(lon0,lat0,lon,lat,rr1) ! GC distance from center 190 CALL dist_lonlat(lon0,lat0,lon_i(n),lat_i(n),rr1) ! GC distance from center 192 191 rr1 = radius*rr1 193 192 IF ( rr1 .LT. R0 ) then … … 212 211 DO i=ii_begin-1,ii_end+1 213 212 n=(j-1)*iim+i 214 CALL xyz2lonlat(xyz_i(n,:),lon,lat) 215 CALL dist_lonlat(lonc1,latc1,lon,lat,rr1) ! GC distance from center 213 CALL dist_lonlat(lonc1,latc1,lon_i(n),lat_i(n),rr1) ! GC distance from center 216 214 rr1 = radius*rr1 217 CALL dist_lonlat(lonc2,latc2,lon ,lat,rr2) ! GC distance from center215 CALL dist_lonlat(lonc2,latc2,lon_i(n),lat_i(n),rr2) ! GC distance from center 218 216 rr2 = radius*rr2 219 217 dd1t1 = rr1/rt … … 244 242 DO i=ii_begin-1,ii_end+1 245 243 n=(j-1)*iim+i 246 CALL xyz2lonlat(xyz_i(n,:),lon,lat) 247 CALL dist_lonlat(lonc1,latc1,lon,lat,rr1) ! GC distance from center 244 CALL dist_lonlat(lonc1,latc1,lon_i(n),lat_i(n),rr1) ! GC distance from center 248 245 rr1 = radius*rr1 249 CALL dist_lonlat(lonc2,latc2,lon ,lat,rr2) ! GC distance from center246 CALL dist_lonlat(lonc2,latc2,lon_i(n),lat_i(n),rr2) ! GC distance from center 250 247 rr2 = radius*rr2 251 248 dd1t1 = rr1/rt -
codes/icosagcm/trunk/src/etat0_dcmip2.f90
r186 r286 86 86 DO i=ii_begin,ii_end 87 87 ij=(j-1)*iim+i 88 CALL comp_all( xyz_i(ij,:), ps(ij),phis(ij),temp(ij,l), dummy1,dummy2)88 CALL comp_all(lon_i(ij),lat_i(ij), ps(ij),phis(ij),temp(ij,l), dummy1,dummy2) 89 89 END DO 90 90 END DO … … 100 100 DO i=ii_begin,ii_end 101 101 ij=(j-1)*iim+i 102 CALL comp_all(xyz_e(ij+u_right,:), dummy1,dummy2,dummy3, ulon(ij+u_right,l),ulat(ij+u_right,l)) 103 CALL comp_all(xyz_e(ij+u_lup,:), dummy1,dummy2,dummy3, ulon(ij+u_lup,l),ulat(ij+u_lup,l)) 104 CALL comp_all(xyz_e(ij+u_ldown,:), dummy1,dummy2,dummy3, ulon(ij+u_ldown,l),ulat(ij+u_ldown,l)) 102 CALL comp_all(lon_e(ij+u_right), lat_e(ij+u_right), & 103 dummy1,dummy2,dummy3, ulon(ij+u_right,l),ulat(ij+u_right,l)) 104 CALL comp_all(lon_e(ij+u_lup), lat_e(ij+u_lup), & 105 dummy1,dummy2,dummy3, ulon(ij+u_lup,l),ulat(ij+u_lup,l)) 106 CALL comp_all(lon_e(ij+u_ldown), lat_e(ij+u_ldown), & 107 dummy1,dummy2,dummy3, ulon(ij+u_ldown,l),ulat(ij+u_ldown,l)) 105 108 END DO 106 109 END DO … … 112 115 CONTAINS 113 116 114 SUBROUTINE comp_all( xyz, psj,phisj,tempj, ulonj,ulatj)117 SUBROUTINE comp_all(lon,lat, psj,phisj,tempj, ulonj,ulatj) 115 118 USE dcmip_initial_conditions_test_1_2_3 116 REAL(rstd), INTENT(IN) :: xyz(3)119 REAL(rstd), INTENT(IN) :: lon, lat 117 120 REAL(rstd), INTENT(OUT) :: psj,phisj,tempj,ulonj,ulatj 118 REAL :: lon,lat,dummy121 REAL :: dummy 119 122 dummy=0. 120 CALL xyz2lonlat(xyz,lon,lat)121 123 SELECT CASE (icase) 122 124 CASE(0) -
codes/icosagcm/trunk/src/etat0_dcmip3.f90
r186 r286 83 83 REAL(rstd) :: Rd ! gas constant of dry air, P=rho.Rd.T 84 84 REAL(rstd) :: alpha, rm 85 REAL(rstd) :: lon,lat,C0, C1, GG85 REAL(rstd) :: C0, C1, GG 86 86 REAL(rstd) :: p0psk, pspsk,r,zz, thetab, thetap 87 87 REAL(rstd) :: dummy, pp … … 103 103 DO i=ii_begin,ii_end 104 104 ij=(j-1)*iim+i 105 CALL xyz2lonlat(xyz_i(ij,:),lon,lat)106 105 107 106 IF(use_dcmip_routine) THEN 108 CALL test3_gravity_wave(lon ,lat,pp,dummy,0, dummy,dummy,dummy,dummy,phis(ij),ps(ij),dummy,dummy)107 CALL test3_gravity_wave(lon_i(ij),lat_i(ij),pp,dummy,0, dummy,dummy,dummy,dummy,phis(ij),ps(ij),dummy,dummy) 109 108 ELSE 110 C1=C0*(cos(2*lat )-1)109 C1=C0*(cos(2*lat_i(ij))-1) 111 110 112 111 !--- GROUND GEOPOTENTIAL … … 120 119 121 120 122 r=radius*acos(sin(latc)*sin(lat )+cos(latc)*cos(lat)*cos(lon-lonc))121 r=radius*acos(sin(latc)*sin(lat_i(ij))+cos(latc)*cos(lat_i(ij))*cos(lon_i(ij)-lonc)) 123 122 s(ij)= d**2/(d**2+r**2) 124 123 END IF … … 134 133 pp=0.5*(p(ij,l+1)+p(ij,l)) ! full-layer pressure 135 134 IF(use_dcmip_routine) THEN 136 CALL xyz2lonlat(xyz_i(ij,:),lon,lat)137 CALL test3_gravity_wave(lon,lat,pp,dummy,0,dummy,dummy,dummy,T(ij,l),dummy,dummy,dummy,dummy)135 CALL test3_gravity_wave(lon_i(ij),lat_i(ij),pp,dummy,0, & 136 dummy,dummy,dummy,T(ij,l),dummy,dummy,dummy,dummy) 138 137 ELSE 139 138 pspsk=(pp/ps(ij))**kappa … … 162 161 ij=(j-1)*iim+i 163 162 IF(use_dcmip_routine) THEN 164 CALL xyz2lonlat(xyz_e(ij+u_right,:),lon,lat) 165 CALL test3_gravity_wave(lon,lat,pp,0.,0, & 166 ulon(ij+u_right,l),ulat(ij+u_right,l),& 163 CALL test3_gravity_wave(lon_e(ij+u_right),lat_e(ij+u_right), & 164 pp,0.,0, ulon(ij+u_right,l),ulat(ij+u_right,l),& 167 165 dummy,dummy,dummy,dummy,dummy,dummy) 168 CALL xyz2lonlat(xyz_e(ij+u_lup,:),lon,lat) 169 CALL test3_gravity_wave(lon,lat,pp,0.,0,& 170 ulon(ij+u_lup,l),ulat(ij+u_lup,l),& 166 CALL test3_gravity_wave(lon_e(ij+u_lup),lat_e(ij+u_lup), & 167 pp,0.,0, ulon(ij+u_lup,l),ulat(ij+u_lup,l),& 171 168 dummy,dummy,dummy,dummy,dummy,dummy) 172 CALL xyz2lonlat(xyz_e(ij+u_ldown,:),lon,lat) 173 CALL test3_gravity_wave(lon,lat,pp,0.,0,& 174 ulon(ij+u_ldown,l),ulat(ij+u_ldown,l),& 169 CALL test3_gravity_wave(lon_e(ij+u_ldown),lat_e(ij+u_ldown), & 170 pp,0.,0, ulon(ij+u_ldown,l),ulat(ij+u_ldown,l),& 175 171 dummy,dummy,dummy,dummy,dummy,dummy) 176 172 ELSE 177 CALL xyz2lonlat(xyz_e(ij+u_right,:),lon,lat) 178 ulon(ij+u_right,l)=u0*cos(lat) 179 ulat(ij+u_right,l)=0 180 181 CALL xyz2lonlat(xyz_e(ij+u_lup,:),lon,lat) 182 ulon(ij+u_lup,l)=u0*cos(lat) 183 ulat(ij+u_lup,l)=0 184 185 CALL xyz2lonlat(xyz_e(ij+u_ldown,:),lon,lat) 186 ulon(ij+u_ldown,l)=u0*cos(lat) 187 ulat(ij+u_ldown,l)=0 173 ulon(ij+u_right,l) = u0*cos(lat_e(ij+u_right)) 174 ulat(ij+u_right,l) = 0 175 ulon(ij+u_lup,l) = u0*cos(lat_e(ij+u_lup)) 176 ulat(ij+u_lup,l) = 0 177 ulon(ij+u_ldown,l) = u0*cos(lat_e(ij+u_ldown)) 178 ulat(ij+u_ldown,l) = 0 188 179 END IF 189 180 ENDDO -
codes/icosagcm/trunk/src/etat0_dcmip4.f90
r186 r286 121 121 ij=(j-1)*iim+i 122 122 123 CALL xyz2lonlat(xyz_e(ij+u_right,:)/radius,lon,lat)123 lon=lon_e(ij+u_right) ; lat=lat_e(ij+u_right) 124 124 K=sin(latc)*sin(lat)+cos(latc)*cos(lat)*cos(lon-lonc) 125 125 r=radius*acos(K) … … 127 127 u(ij+u_right,l) = utot * sum(elon_e(ij+u_right,:) * ep_e(ij+u_right,:)) 128 128 129 130 CALL xyz2lonlat(xyz_e(ij+u_lup,:)/radius,lon,lat) 129 lon=lon_e(ij+u_lup) ; lat=lat_e(ij+u_lup) 131 130 K=sin(latc)*sin(lat)+cos(latc)*cos(lat)*cos(lon-lonc) 132 131 r=radius*acos(K) … … 134 133 u(ij+u_lup,l) = utot * sum(elon_e(ij+u_lup,:) * ep_e(ij+u_lup,:)) 135 134 136 CALL xyz2lonlat(xyz_e(ij+u_ldown,:)/radius,lon,lat)135 lon=lon_e(ij+u_ldown) ; lat=lat_e(ij+u_ldown) 137 136 K=sin(latc)*sin(lat)+cos(latc)*cos(lat)*cos(lon-lonc) 138 137 r=radius*acos(K) … … 151 150 DO i=ii_begin-1,ii_end+1 152 151 ij=(j-1)*iim+i 153 CALL xyz2lonlat(xyz_i(ij,:)/radius,lon,lat) 154 152 lat=lat_i(ij) 155 153 Y(ij,l)=((-2*sin(lat)**6*(cos(lat)**2+1./3)+10./63)*2*u0*cos(etav(l))**1.5 & 156 154 + (8./5*cos(lat)**3*(sin(lat)**2+2./3)-Pi/4)*radius*Omega) … … 168 166 DO i=ii_begin,ii_end 169 167 ij=(j-1)*iim+i 170 CALL xyz2lonlat(xyz_i(ij,:)/radius,lon,lat)168 lat=lat_i(ij) 171 169 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 & 172 170 +(8./5*cos(lat)**3 * (sin(lat)**2 + 2./3) - Pi/4)*radius*Omega ) … … 191 189 DO i=ii_begin,ii_end 192 190 ij=(j-1)*iim+i 193 CALL xyz2lonlat(xyz_i(ij,:)/radius,lon,lat)191 lon=lon_i(ij) ; lat=lat_i(ij) 194 192 dthetaodeta=dthetaodeta_ave + 3./4. * Pi * u0/Rd*(1-kappa)*eta(l)**(-kappa)*sin(etav(l))*cos(etav(l))**0.5 * Y(ij,l) & 195 193 + 3/8. * Pi**2*u0/Rd * eta(l)**(1-kappa) * cos(etav(l))**1.5 * Y(ij,l) & … … 220 218 DO i=ii_begin,ii_end 221 219 ij=(j-1)*iim+i 222 CALL xyz2lonlat(xyz_i(ij,:)/radius,lon,lat) 223 q(ij,l,1)=q0*exp(-(lat/latw)**4)*exp(-((eta(l)-1)*preff/pw)**2) 220 q(ij,l,1)=q0*exp(-(lat_i(ij)/latw)**4)*exp(-((eta(l)-1)*preff/pw)**2) 224 221 ENDDO 225 222 ENDDO -
codes/icosagcm/trunk/src/etat0_heldsz.f90
r186 r286 6 6 TYPE(t_field),POINTER :: f_theta_eq(:) 7 7 TYPE(t_field),POINTER :: f_theta(:) 8 TYPE(t_field),POINTER :: f_clat(:) ! FIXME, duplication9 8 10 9 REAL(rstd),ALLOCATABLE,SAVE :: knewt_t(:),kfrict(:) … … 65 64 REAL(rstd),POINTER :: u(:,:) 66 65 REAL(rstd),POINTER :: q(:,:,:) 67 REAL(rstd),POINTER :: clat(:)68 66 REAL(rstd),POINTER :: theta_eq(:,:) 69 67 REAL(rstd),POINTER :: theta(:,:) … … 78 76 79 77 theta_eq=f_theta_eq(ind) 80 clat=f_clat(ind) 81 CALL compute_Teq(clat,theta_eq) ! FIXME : already done by Init_Teq 78 CALL compute_Teq(lat_i,theta_eq) ! FIXME : already done by Init_Teq 82 79 83 80 ps=f_ps(ind) … … 111 108 CALL allocate_field(f_theta,field_t,type_real,llm) 112 109 CALL allocate_field(f_theta_eq,field_t,type_real,llm) 113 CALL allocate_field(f_clat,field_t,type_real)114 110 ALLOCATE(knewt_t(llm)); ALLOCATE( kfrict(llm)) 115 111 … … 147 143 CALL swap_dimensions(ind) 148 144 CALL swap_geometry(ind) 149 clat=f_clat(ind)150 145 theta_eq=f_theta_eq(ind) 151 CALL compute_Teq( clat,theta_eq)146 CALL compute_Teq(lat_i,theta_eq) 152 147 ENDDO 153 148 … … 159 154 END SUBROUTINE init_Teq 160 155 161 SUBROUTINE compute_Teq( clat,theta_eq)156 SUBROUTINE compute_Teq(lat,theta_eq) 162 157 USE icosa 163 158 USE disvert_mod 164 159 IMPLICIT NONE 165 REAL(rstd),INTENT( OUT) :: clat(iim*jjm)160 REAL(rstd),INTENT(IN) :: lat(iim*jjm) 166 161 REAL(rstd),INTENT(OUT) :: theta_eq(iim*jjm,llm) 167 162 168 REAL(rstd) :: lon, lat, r, zsig, ddsin, tetastrat, tetajl 169 REAL(rstd) :: slat(iim*jjm) 163 REAL(rstd) :: r, zsig, ddsin, tetastrat, tetajl 170 164 INTEGER :: i,j,l,ij 171 172 DO j=jj_begin-1,jj_end+1173 DO i=ii_begin-1,ii_end+1174 ij=(j-1)*iim+i175 CALL xyz2lonlat(xyz_i(ij,:),lon,lat)176 clat(ij)=cos(lat)177 slat(ij)=sin(lat)178 ENDDO179 ENDDO180 165 181 166 DO l=1,llm … … 185 170 DO i=ii_begin-1,ii_end+1 186 171 ij=(j-1)*iim+i 187 ddsin= slat(ij)172 ddsin=SIN(lat(ij)) 188 173 tetajl=teta0-delt_y*ddsin*ddsin+eps*ddsin & 189 174 -delt_z*(1.-ddsin*ddsin)*log(zsig) … … 244 229 theta_eq=f_theta_eq(ind) 245 230 theta=f_theta(ind) 246 clat=f_clat(ind) 247 CALL compute_heldsz(ps,theta_eq,clat, theta_rhodz,u, theta) 231 CALL compute_heldsz(ps,theta_eq,lat_i, theta_rhodz,u, theta) 248 232 ENDDO 249 233 END SUBROUTINE held_suarez 250 234 251 SUBROUTINE compute_heldsz(ps,theta_eq, clat, theta_rhodz,u, theta)235 SUBROUTINE compute_heldsz(ps,theta_eq,lat, theta_rhodz,u, theta) 252 236 USE icosa 253 237 USE theta2theta_rhodz_mod … … 255 239 REAL(rstd),INTENT(IN) :: ps(iim*jjm) 256 240 REAL(rstd),INTENT(IN) :: theta_eq(iim*jjm,llm) 257 REAL(rstd),INTENT(IN) :: clat(iim*jjm)241 REAL(rstd),INTENT(IN) :: lat(iim*jjm) 258 242 REAL(rstd),INTENT(INOUT) :: theta_rhodz(iim*jjm,llm) 259 243 REAL(rstd),INTENT(INOUT) :: u(3*iim*jjm,llm) … … 268 252 ij=(j-1)*iim+i 269 253 theta(ij,l)=theta(ij,l) - dt*(theta(ij,l)-theta_eq(ij,l))* & 270 (knewt_g+knewt_t(l)* clat(ij)**4 )254 (knewt_g+knewt_t(l)*COS(lat(ij))**4 ) 271 255 ENDDO 272 256 ENDDO -
codes/icosagcm/trunk/src/geometry.f90
r266 r286 3 3 4 4 TYPE t_geometry 5 TYPE(t_field),POINTER :: centroid(:) 5 6 TYPE(t_field),POINTER :: xyz_i(:) 6 TYPE(t_field),POINTER :: centroid(:)7 7 TYPE(t_field),POINTER :: xyz_e(:) 8 8 TYPE(t_field),POINTER :: xyz_v(:) 9 TYPE(t_field),POINTER :: lon_i(:) 10 TYPE(t_field),POINTER :: lon_e(:) 11 TYPE(t_field),POINTER :: lat_i(:) 12 TYPE(t_field),POINTER :: lat_e(:) 9 13 TYPE(t_field),POINTER :: ep_e(:) 10 14 TYPE(t_field),POINTER :: et_e(:) … … 30 34 REAL(rstd),POINTER :: Ai(:) ! area of a cell 31 35 !$OMP THREADPRIVATE(Ai) 36 REAL(rstd),POINTER :: centroid(:,:) ! coordinate of the centroid of the cell 37 !$OMP THREADPRIVATE(centroid) 32 38 REAL(rstd),POINTER :: xyz_i(:,:) ! coordinate of the center of the cell (voronoi) 33 39 !$OMP THREADPRIVATE(xyz_i) 34 REAL(rstd),POINTER :: centroid(:,:) ! coordinate of the centroid of the cell35 !$OMP THREADPRIVATE(centroid)36 40 REAL(rstd),POINTER :: xyz_e(:,:) ! coordinate of a wind point on the cell on a edge 37 41 !$OMP THREADPRIVATE(xyz_e) 42 REAL(rstd),POINTER :: xyz_v(:,:) ! coordinate of a vertex (center of the dual mesh) 43 !$OMP THREADPRIVATE(xyz_v) 44 REAL(rstd),POINTER :: lon_i(:) ! longitude of the center of the cell (voronoi) 45 !$OMP THREADPRIVATE(lon_i) 46 REAL(rstd),POINTER :: lon_e(:) ! longitude of a wind point on the cell on a edge 47 !$OMP THREADPRIVATE(lon_e) 48 REAL(rstd),POINTER :: lat_i(:) ! latitude of the center of the cell (voronoi) 49 !$OMP THREADPRIVATE(lat_i) 50 REAL(rstd),POINTER :: lat_e(:) ! latitude of a wind point on the cell on a edge 51 !$OMP THREADPRIVATE(lat_e) 38 52 REAL(rstd),POINTER :: ep_e(:,:) ! perpendicular unit vector of a edge (outsider) 39 53 !$OMP THREADPRIVATE(ep_e) … … 48 62 REAL(rstd),POINTER :: elat_e(:,:) ! unit latitude vector on a wind point 49 63 !$OMP THREADPRIVATE(elat_e) 50 REAL(rstd),POINTER :: xyz_v(:,:) ! coordinate of a vertex (center of the dual mesh)51 !$OMP THREADPRIVATE(xyz_v)52 64 REAL(rstd),POINTER :: Av(:) ! area of dual mesk cell 53 65 !$OMP THREADPRIVATE(Av) … … 84 96 85 97 CALL allocate_field(geom%Ai,field_t,type_real,name='Ai') 98 86 99 CALL allocate_field(geom%xyz_i,field_t,type_real,3) 100 CALL allocate_field(geom%lon_i,field_t,type_real) 101 CALL allocate_field(geom%lat_i,field_t,type_real) 102 CALL allocate_field(geom%elon_i,field_t,type_real,3) 103 CALL allocate_field(geom%elat_i,field_t,type_real,3) 87 104 CALL allocate_field(geom%centroid,field_t,type_real,3) 105 88 106 CALL allocate_field(geom%xyz_e,field_u,type_real,3) 107 CALL allocate_field(geom%lon_e,field_u,type_real) 108 CALL allocate_field(geom%lat_e,field_u,type_real) 109 CALL allocate_field(geom%elon_e,field_u,type_real,3) 110 CALL allocate_field(geom%elat_e,field_u,type_real,3) 89 111 CALL allocate_field(geom%ep_e,field_u,type_real,3) 90 112 CALL allocate_field(geom%et_e,field_u,type_real,3) 91 CALL allocate_field(geom%elon_i,field_t,type_real,3) 92 CALL allocate_field(geom%elat_i,field_t,type_real,3) 93 CALL allocate_field(geom%elon_e,field_u,type_real,3) 94 CALL allocate_field(geom%elat_e,field_u,type_real,3) 113 95 114 CALL allocate_field(geom%xyz_v,field_z,type_real,3) 96 115 CALL allocate_field(geom%de,field_u,type_real) … … 119 138 ep_e=geom%ep_e(ind) 120 139 et_e=geom%et_e(ind) 140 lon_i=geom%lon_i(ind) 141 lat_i=geom%lat_i(ind) 142 lon_e=geom%lon_e(ind) 143 lat_e=geom%lat_e(ind) 121 144 elon_i=geom%elon_i(ind) 122 145 elat_i=geom%elat_i(ind) … … 394 417 395 418 CALL xyz2lonlat(xyz_i(n,:),lon,lat) 419 lon_i(n)=lon 420 lat_i(n)=lat 396 421 elon_i(n,1) = -sin(lon) 397 422 elon_i(n,2) = cos(lon) … … 403 428 404 429 CALL xyz2lonlat(xyz_e(n+u_right,:),lon,lat) 430 lon_e(n+u_right)=lon 431 lat_e(n+u_right)=lat 405 432 elon_e(n+u_right,1) = -sin(lon) 406 433 elon_e(n+u_right,2) = cos(lon) … … 411 438 412 439 CALL xyz2lonlat(xyz_e(n+u_lup,:),lon,lat) 440 lon_e(n+u_lup)=lon 441 lat_e(n+u_lup)=lat 413 442 elon_e(n+u_lup,1) = -sin(lon) 414 443 elon_e(n+u_lup,2) = cos(lon) … … 419 448 420 449 CALL xyz2lonlat(xyz_e(n+u_ldown,:),lon,lat) 450 lon_e(n+u_ldown)=lon 451 lat_e(n+u_ldown)=lat 421 452 elon_e(n+u_ldown,1) = -sin(lon) 422 453 elon_e(n+u_ldown,2) = cos(lon) -
codes/icosagcm/trunk/src/physics_interface.f90
r214 r286 194 194 CALL swap_geometry(ind) 195 195 CALL pack_domain_2D(pack_info(ind), Ai, physics_inout%Ai) 196 CALL pack_lonlat(pack_info(ind)) 196 CALL pack_domain_2D(pack_info(ind), lon_i, physics_inout%lon) 197 CALL pack_domain_2D(pack_info(ind), lat_i, physics_inout%lat) 197 198 END DO 198 199 END SUBROUTINE init_pack_after 199 200 SUBROUTINE pack_lonlat(info)201 USE icosa202 IMPLICIT NONE203 TYPE(t_pack_info) :: info204 REAL(rstd) :: lon(iim*jjm), lat(iim*jjm)205 INTEGER :: ij206 DO ij=1,iim*jjm207 CALL xyz2lonlat(xyz_i(ij,:),lon(ij),lat(ij))208 ENDDO209 CALL pack_domain_2D(info, lon, physics_inout%lon)210 CALL pack_domain_2D(info, lat, physics_inout%lat)211 END SUBROUTINE pack_lonlat212 200 213 201 !-------------------------------- Pack / Unpack 2D ---------------------------
Note: See TracChangeset
for help on using the changeset viewer.