Changeset 985 for LMDZ4/trunk/libf/dyn3dpar/filtreg_p.F
- Timestamp:
- Jul 30, 2008, 5:50:03 PM (16 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
LMDZ4/trunk/libf/dyn3dpar/filtreg_p.F
r763 r985 1 2 1 3 SUBROUTINE filtreg_p ( champ, ibeg, iend, nlat, nbniv, 2 4 . ifiltre, iaire, griscal ,iter) 3 5 USE Parallel, only : OMP_CHUNK 6 USE mod_filtre_fft 7 USE timer_filtre 4 8 IMPLICIT NONE 5 9 … … 59 63 , , matricevn(iim,iim,nfilvn),matricevs(iim,iim,nfilvs) 60 64 , , matrinvn(iim,iim,nfilun),matrinvs (iim,iim,nfilus) 61 REAL eignq(iim), sdd1(iim),sdd2(iim) 65 cym REAL eignq(iim), sdd1(iim),sdd2(iim) 66 67 REAL eignq(iim) 68 REAL :: sdd1(iim),sdd2(iim) 69 62 70 LOGICAL griscal 63 71 INTEGER hemisph, iaire 64 c 72 73 REAL :: champ_fft(iip1,nlat,nbniv) 74 REAL :: champ_in(iip1,nlat,nbniv) 75 76 REAL,SAVE,TARGET :: sddu_loc(iim) 77 REAL,SAVE,TARGET :: sddv_loc(iim) 78 REAL,SAVE,TARGET :: unsddu_loc(iim) 79 REAL,SAVE,TARGET :: unsddv_loc(iim) 80 c$OMP THREADPRIVATE(sddu_loc,sddv_loc,unsddu_loc,unsddv_loc) 81 LOGICAL,SAVE :: first=.TRUE. 82 c$OMP THREADPRIVATE(first) 83 84 IF (first) THEN 85 sddu_loc(1:iim)=sddu(1:iim) 86 sddv_loc(1:iim)=sddv(1:iim) 87 unsddu_loc(1:iim)=unsddu(1:iim) 88 unsddv_loc(1:iim)=unsddv(1:iim) 89 CALL Init_timer 90 first=.FALSE. 91 c PRINT *,"----> sddu_loc=",sddu_loc 92 c PRINT *,"----> sddv_loc=",sddv_loc 93 c PRINT *,"----> unsddu_loc=",unsddu_loc 94 c PRINT *,"----> unsddv_loc=",unsddv_loc 95 ENDIF 96 97 c$OMP MASTER 98 CALL start_timer 99 c$OMP END MASTER 65 100 66 101 IF(ifiltre.EQ.1.or.ifiltre.EQ.-1) … … 97 132 c 98 133 IF( iaire.EQ.1 ) THEN 99 CALL SCOPY( iim, sddv, 1, sdd1, 1 ) 100 CALL SCOPY( iim, unsddv, 1, sdd2, 1 ) 134 cym CALL SCOPY( iim, sddv, 1, sdd1, 1 ) 135 cym CALL SCOPY( iim, unsddv, 1, sdd2, 1 ) 136 cym sdd1=>sddv_loc 137 cym sdd2=>unsddv_loc 138 sdd1(1:iim)=sddv_loc(1:iim) 139 sdd2(1:iim)=unsddv_loc(1:iim) 101 140 ELSE 102 CALL SCOPY( iim, unsddv, 1, sdd1, 1 ) 103 CALL SCOPY( iim, sddv, 1, sdd2, 1 ) 141 cym CALL SCOPY( iim, unsddv, 1, sdd1, 1 ) 142 cym CALL SCOPY( iim, sddv, 1, sdd2, 1 ) 143 sdd1(1:iim)=unsddv_loc(1:iim) 144 sdd2(1:iim)=sddv_loc(1:iim) 104 145 END IF 105 146 c … … 116 157 c 117 158 IF( iaire.EQ.1 ) THEN 118 CALL SCOPY( iim, sddu, 1, sdd1, 1 ) 119 CALL SCOPY( iim, unsddu, 1, sdd2, 1 ) 159 cym CALL SCOPY( iim, sddu, 1, sdd1, 1 ) 160 cym CALL SCOPY( iim, unsddu, 1, sdd2, 1 ) 161 cym sdd1=>sddu_loc 162 cym sdd2=>unsddu_loc 163 sdd1(1:iim)=sddu_loc(1:iim) 164 sdd2(1:iim)=unsddu_loc(1:iim) 165 120 166 ELSE 121 CALL SCOPY( iim, unsddu, 1, sdd1, 1 ) 122 CALL SCOPY( iim, sddu, 1, sdd2, 1 ) 167 cym CALL SCOPY( iim, unsddu, 1, sdd1, 1 ) 168 cym CALL SCOPY( iim, sddu, 1, sdd2, 1 ) 169 cym sdd1=>unsddu_loc 170 cym sdd2=>sddu_loc 171 sdd1(1:iim)=unsddu_loc(1:iim) 172 sdd2(1:iim)=sddu_loc(1:iim) 123 173 END IF 124 174 c … … 129 179 END IF 130 180 END IF 181 182 c PRINT *,"APPEL a filtreg --> sdd1=",sdd1 183 c PRINT *,"APPEL a filtreg --> sdd2=",sdd2 184 c PRINT *,"----> sddu_loc=",sddu_loc 185 c PRINT *,"----> sddv_loc=",sddv_loc 186 c PRINT *,"----> unsddu_loc=",unsddu_loc 187 c PRINT *,"----> unsddv_loc=",unsddv_loc 188 131 189 c 132 190 c … … 143 201 END IF 144 202 203 204 cccccccccccccccccccccccccccccccccccccccccccc 205 c Utilisation du filtre classique 206 cccccccccccccccccccccccccccccccccccccccccccc 207 208 IF (.NOT. use_filtre_fft) THEN 209 145 210 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 146 211 DO 50 l = 1, nbniv 147 DO 30 j = jdfil,jffil212 DO 30 j = jdfil,jffil 148 213 149 214 150 DO 5 i = 1, iim 151 champ(i,j,l) = champ(i,j,l) * sdd1(i) 152 5 CONTINUE 153 c 154 155 IF( hemisph. EQ. 1 ) THEN 156 157 IF( ifiltre. EQ. -2 ) THEN 158 #ifdef CRAY 159 CALL MXVA( matrinvn(1,1,j), 1, iim, champ(1,j,l), 1, eignq , 160 * 1, iim, iim ) 161 #else 162 #ifdef BLAS 163 CALL SGEMV("N", iim,iim, 1.0, matrinvn(1,1,j),iim, 164 . champ(1,j,l), 1, 0.0, eignq, 1) 165 #else 166 DO k = 1, iim 167 eignq(k) = 0.0 168 ENDDO 169 DO k = 1, iim 170 DO i = 1, iim 171 eignq(k) = eignq(k) + matrinvn(k,i,j)*champ(i,j,l) 172 ENDDO 173 ENDDO 174 #endif 175 #endif 176 ELSE IF ( griscal ) THEN 177 #ifdef CRAY 178 CALL MXVA( matriceun(1,1,j), 1, iim, champ(1,j,l), 1, eignq , 179 * 1, iim, iim ) 180 #else 181 #ifdef BLAS 182 CALL SGEMV("N", iim,iim, 1.0, matriceun(1,1,j),iim, 183 . champ(1,j,l), 1, 0.0, eignq, 1) 184 #else 185 DO k = 1, iim 186 eignq(k) = 0.0 187 ENDDO 188 DO i = 1, iim 189 DO k = 1, iim 190 eignq(k) = eignq(k) + matriceun(k,i,j)*champ(i,j,l) 191 ENDDO 192 ENDDO 193 #endif 194 #endif 195 ELSE 196 #ifdef CRAY 197 CALL MXVA( matricevn(1,1,j), 1, iim, champ(1,j,l), 1, eignq , 198 * 1, iim, iim ) 199 #else 200 #ifdef BLAS 201 CALL SGEMV("N", iim,iim, 1.0, matricevn(1,1,j),iim, 202 . champ(1,j,l), 1, 0.0, eignq, 1) 203 #else 204 DO k = 1, iim 205 eignq(k) = 0.0 206 ENDDO 207 DO i = 1, iim 208 DO k = 1, iim 209 eignq(k) = eignq(k) + matricevn(k,i,j)*champ(i,j,l) 210 ENDDO 211 ENDDO 212 #endif 213 #endif 214 ENDIF 215 216 ELSE 217 218 IF( ifiltre. EQ. -2 ) THEN 219 #ifdef CRAY 220 CALL MXVA( matrinvs(1,1,j-jfiltsu+1), 1, iim, champ(1,j,l),1 , 221 * eignq, 1, iim, iim ) 222 #else 223 #ifdef BLAS 224 CALL SGEMV("N", iim,iim, 1.0, matrinvs(1,1,j-jfiltsu+1),iim, 225 . champ(1,j,l), 1, 0.0, eignq, 1) 226 #else 227 DO k = 1, iim 228 eignq(k) = 0.0 229 ENDDO 230 DO i = 1, iim 231 DO k = 1, iim 232 eignq(k) = eignq(k) + matrinvs(k,i,j-jfiltsu+1)*champ(i,j,l) 233 ENDDO 234 ENDDO 235 #endif 236 #endif 237 ELSE IF ( griscal ) THEN 238 #ifdef CRAY 239 CALL MXVA( matriceus(1,1,j-jfiltsu+1), 1, iim, champ(1,j,l),1 , 240 * eignq, 1, iim, iim ) 241 #else 242 #ifdef BLAS 243 CALL SGEMV("N", iim,iim, 1.0, matriceus(1,1,j-jfiltsu+1),iim, 244 . champ(1,j,l), 1, 0.0, eignq, 1) 245 #else 246 DO k = 1, iim 247 eignq(k) = 0.0 248 ENDDO 249 DO i = 1, iim 250 DO k = 1, iim 251 eignq(k) = eignq(k) + matriceus(k,i,j-jfiltsu+1)*champ(i,j,l) 252 ENDDO 253 ENDDO 254 #endif 255 #endif 256 ELSE 257 #ifdef CRAY 258 CALL MXVA( matricevs(1,1,j-jfiltsv+1), 1, iim, champ(1,j,l),1 , 259 * eignq, 1, iim, iim ) 260 #else 261 #ifdef BLAS 262 CALL SGEMV("N", iim,iim, 1.0, matricevs(1,1,j-jfiltsv+1),iim, 263 . champ(1,j,l), 1, 0.0, eignq, 1) 264 #else 265 DO k = 1, iim 266 eignq(k) = 0.0 267 ENDDO 268 DO i = 1, iim 269 DO k = 1, iim 270 eignq(k) = eignq(k) + matricevs(k,i,j-jfiltsv+1)*champ(i,j,l) 271 ENDDO 272 ENDDO 273 #endif 274 #endif 275 ENDIF 276 277 ENDIF 278 c 279 IF( ifiltre.EQ. 2 ) THEN 280 DO 15 i = 1, iim 281 champ( i,j,l ) = ( champ(i,j,l) + eignq(i) ) * sdd2(i) 282 15 CONTINUE 283 ELSE 284 DO 16 i=1,iim 285 champ( i,j,l ) = ( champ(i,j,l) - eignq(i) ) * sdd2(i) 286 16 CONTINUE 287 ENDIF 288 c 289 champ( iip1,j,l ) = champ( 1,j,l ) 290 c 291 30 CONTINUE 215 DO 5 i = 1, iim 216 champ(i,j,l) = champ(i,j,l) * sdd1(i) 217 5 CONTINUE 218 c 219 220 IF( hemisph. EQ. 1 ) THEN 221 222 IF( ifiltre. EQ. -2 ) THEN 223 224 225 CALL SGEMV("N", iim,iim, 1.0, matrinvn(1,1,j),iim, 226 . champ(1,j,l), 1, 0.0, eignq, 1) 227 228 229 ELSE IF ( griscal ) THEN 230 231 CALL SGEMV("N", iim,iim, 1.0, matriceun(1,1,j),iim, 232 . champ(1,j,l), 1, 0.0, eignq, 1) 233 234 ELSE 235 236 CALL SGEMV("N", iim,iim, 1.0, matricevn(1,1,j),iim, 237 . champ(1,j,l), 1, 0.0, eignq, 1) 238 ENDIF 239 240 ELSE 241 242 IF( ifiltre. EQ. -2 ) THEN 243 244 CALL SGEMV("N",iim,iim,1.0, matrinvs(1,1,j-jfiltsu+1),iim, 245 . champ(1,j,l), 1, 0.0, eignq, 1) 246 247 ELSE IF ( griscal ) THEN 248 249 CALL SGEMV("N",iim,iim,1.0,matriceus(1,1,j-jfiltsu+1),iim, 250 . champ(1,j,l), 1, 0.0, eignq, 1) 251 ELSE 252 253 CALL SGEMV("N",iim,iim,1.0,matricevs(1,1,j-jfiltsv+1),iim, 254 . champ(1,j,l), 1, 0.0, eignq, 1) 255 ENDIF 256 257 ENDIF 258 259 260 c 261 IF( ifiltre.EQ. 2 ) THEN 262 263 DO 15 i = 1, iim 264 champ( i,j,l ) = ( champ(i,j,l) + eignq(i) ) * sdd2(i) 265 15 CONTINUE 266 267 ELSE 268 269 DO 16 i=1,iim 270 champ( i,j,l ) = ( champ(i,j,l) - eignq(i) ) * sdd2(i) 271 16 CONTINUE 272 273 ENDIF 274 c 275 champ( iip1,j,l ) = champ( 1,j,l ) 276 c 277 30 CONTINUE 292 278 c 293 279 50 CONTINUE 294 280 c$OMP END DO NOWAIT 295 c 281 282 ccccccccccccccccccccccccccccccccccccccccccccc 283 c Utilisation du filtre FFT 284 ccccccccccccccccccccccccccccccccccccccccccccc 285 286 ELSE 287 288 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 289 DO l=1,nbniv 290 DO j=jdfil,jffil 291 DO i = 1, iim 292 champ( i,j,l)= champ(i,j,l)*sdd1(i) 293 champ_fft( i,j,l) = champ(i,j,l) 294 ENDDO 295 ENDDO 296 ENDDO 297 c$OMP END DO NOWAIT 298 299 IF (jdfil<=jffil) THEN 300 IF( ifiltre. EQ. -2 ) THEN 301 CALL Filtre_inv_fft(champ_fft,nlat,jdfil,jffil,nbniv) 302 ELSE IF ( griscal ) THEN 303 CALL Filtre_u_fft(champ_fft,nlat,jdfil,jffil,nbniv) 304 ELSE 305 CALL Filtre_v_fft(champ_fft,nlat,jdfil,jffil,nbniv) 306 ENDIF 307 ENDIF 308 309 310 IF( ifiltre.EQ. 2 ) THEN 311 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 312 DO l=1,nbniv 313 DO j=jdfil,jffil 314 DO i = 1, iim 315 champ( i,j,l)=(champ(i,j,l)+champ_fft(i,j,l)) 316 . *sdd2(i) 317 ENDDO 318 ENDDO 319 ENDDO 320 c$OMP END DO NOWAIT 321 ELSE 322 323 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 324 DO l=1,nbniv 325 DO j=jdfil,jffil 326 DO i = 1, iim 327 champ(i,j,l)=(champ(i,j,l)-champ_fft(i,j,l)) 328 . *sdd2(i) 329 ENDDO 330 ENDDO 331 ENDDO 332 c$OMP END DO NOWAIT 333 ENDIF 334 c 335 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 336 DO l=1,nbniv 337 DO j=jdfil,jffil 338 ! champ_FFT( iip1,j,l ) = champ_FFT( 1,j,l ) 339 champ( iip1,j,l ) = champ( 1,j,l ) 340 ENDDO 341 ENDDO 342 c$OMP END DO NOWAIT 343 ENDIF 344 c Fin de la zone de filtrage 345 346 296 347 100 CONTINUE 348 349 ! DO j=1,nlat 350 ! 351 ! PRINT *,"check FFT ----> Delta(",j,")=", 352 ! & sum(champ(:,j,:)-champ_fft(:,j,:))/sum(champ(:,j,:)), 353 ! & sum(champ(:,j,:)-champ_in(:,j,:))/sum(champ(:,j,:)) 354 ! ENDDO 355 356 ! PRINT *,"check FFT ----> Delta(",j,")=", 357 ! & sum(champ-champ_fft)/sum(champ) 358 ! 359 297 360 c 298 361 1111 FORMAT(//20x,'ERREUR dans le dimensionnement du tableau CHAMP a … … 300 363 2222 FORMAT(//20x,'ERREUR dans le dimensionnement du tableau CHAMP a fi 301 364 *ltrer, sur la grille de V ou de Z'/) 365 c$OMP MASTER 366 CALL stop_timer 367 c$OMP END MASTER 302 368 RETURN 303 369 END
Note: See TracChangeset
for help on using the changeset viewer.