Changeset 297 for LMDZ.3.3/trunk/libf
- Timestamp:
- Nov 29, 2001, 1:09:18 PM (23 years ago)
- File:
-
- 1 edited
-
LMDZ.3.3/trunk/libf/dyn3d/fxhyp.F (modified) (17 diffs)
Legend:
- Unmodified
- Added
- Removed
-
LMDZ.3.3/trunk/libf/dyn3d/fxhyp.F
r251 r297 1 c 2 c $Header$ 3 c 1 4 SUBROUTINE fxhyp ( xzoomdeg,grossism,dzoom,tau , 2 5 , rlonm025,xprimm025,rlonv,xprimv,rlonu,xprimu,rlonp025,xprimp025, … … 20 23 INTEGER nmax, nmax2 21 24 PARAMETER ( nmax = 30000, nmax2 = 2*nmax ) 25 c 26 LOGICAL scal180 27 PARAMETER ( scal180 = .TRUE. ) 28 29 c scal180 = .TRUE. si on veut avoir le premier point scalaire pour 30 c une grille reguliere ( grossism = 1.,tau=0.,clon=0. ) a -180. degres. 31 c sinon scal180 = .FALSE. 22 32 23 33 #include "dimensions.h" … … 42 52 REAL*8 pi,depi,epsilon,xzoom,fa,fb 43 53 REAL*8 Xf1, Xfi , a0,a1,a2,a3,xi2 44 INTEGER i,it,ik,iter,ii,idif 54 INTEGER i,it,ik,iter,ii,idif,ii1,ii2 45 55 REAL*8 xi,xo1,xmoy,xlon2,fxm,Xprimin 46 REAL*8 champmin,champmax 56 REAL*8 champmin,champmax,decalx 47 57 INTEGER is2 48 58 SAVE is2 … … 55 65 epsilon = 1.e-3 56 66 xzoom = xzoomdeg * pi/180. 57 67 c 68 decalx = .75 69 IF( grossism.EQ.1..AND.scal180 ) THEN 70 decalx = 1. 71 ENDIF 72 73 WRITE(6,*) 'FXHYP scal180,decalx', scal180,decalx 74 c 58 75 IF( dzoom.LT.1.) THEN 59 76 dzoom = dzoom * depi … … 73 90 ENDDO 74 91 75 DO i = 0, nmax292 DO i = nmax, nmax2 76 93 77 94 fa = tau* ( dzoom/2. - xtild(i) ) 78 95 fb = xtild(i) * ( pi - xtild(i) ) 79 96 80 81 if ( xtild(i) .ne. 0. .and. xtild(i) .ne. pi) then82 97 IF( 200.* fb .LT. - fa ) THEN 83 98 fhyp ( i) = - 1. … … 85 100 fhyp ( i) = 1. 86 101 ELSE 87 fhyp(i) = TANH ( fa/fb ) 102 IF( ABS(fa).LT.1.e-13.AND.ABS(fb).LT.1.e-13) THEN 103 IF( 200.*fb + fa.LT.1.e-10 ) THEN 104 fhyp ( i ) = - 1. 105 ELSEIF( 200.*fb - fa.LT.1.e-10 ) THEN 106 fhyp ( i ) = 1. 107 ENDIF 108 ELSE 109 fhyp ( i ) = TANH ( fa/fb ) 110 ENDIF 88 111 ENDIF 89 endif 90 91 IF ( xtild(i).EQ. 0. ) fhyp(i) = 1. 92 IF ( xtild(i).EQ. pi ) fhyp(i) = -1. 112 113 IF ( xtild(i).EQ. 0. ) fhyp(i) = 1. 114 IF ( xtild(i).EQ. pi ) fhyp(i) = -1. 93 115 94 116 ENDDO … … 100 122 101 123 DO i = nmax +1,nmax2 124 125 xmoy = 0.5 * ( xtild(i-1) + xtild( i ) ) 126 fa = tau* ( dzoom/2. - xmoy ) 127 fb = xmoy * ( pi - xmoy ) 128 129 IF( 200.* fb .LT. - fa ) THEN 130 fxm = - 1. 131 ELSEIF( 200. * fb .LT. fa ) THEN 132 fxm = 1. 133 ELSE 134 IF( ABS(fa).LT.1.e-13.AND.ABS(fb).LT.1.e-13) THEN 135 IF( 200.*fb + fa.LT.1.e-10 ) THEN 136 fxm = - 1. 137 ELSEIF( 200.*fb - fa.LT.1.e-10 ) THEN 138 fxm = 1. 139 ENDIF 140 ELSE 141 fxm = TANH ( fa/fb ) 142 ENDIF 143 ENDIF 144 145 IF ( xmoy.EQ. 0. ) fxm = 1. 146 IF ( xmoy.EQ. pi ) fxm = -1. 147 148 ffdx = ffdx + fxm * ( xtild(i) - xtild(i-1) ) 149 150 ENDDO 151 152 beta = ( grossism * ffdx - pi ) / ( ffdx - pi ) 153 154 IF( 2.*beta - grossism.LE. 0.) THEN 155 WRITE(6,*) ' ** Attention ! La valeur beta calculee dans la rou 156 ,tine fxhyp est mauvaise ! ' 157 WRITE(6,*)'Modifier les valeurs de grossismx ,tau ou dzoomx ', 158 , ' et relancer ! *** ' 159 CALL ABORT 160 ENDIF 161 c 162 c ..... calcul de Xprimt ..... 163 c 164 165 DO i = nmax, nmax2 166 Xprimt(i) = beta + ( grossism - beta ) * fhyp(i) 167 ENDDO 168 c 169 DO i = nmax+1, nmax2 170 Xprimt( nmax2 - i ) = Xprimt( i ) 171 ENDDO 172 c 173 174 c ..... Calcul de Xf ........ 175 176 Xf(0) = - pi 177 178 DO i = nmax +1, nmax2 102 179 103 180 xmoy = 0.5 * ( xtild(i-1) + xtild( i ) ) … … 113 190 ENDIF 114 191 115 IF ( xmoy.EQ. 0. ) fhyp(i) = 1.116 IF ( xmoy.EQ. pi ) fhyp(i) = -1.117 ffdx = ffdx + fxm * ( xtild(i) - xtild(i-1) )118 119 ENDDO120 121 beta = ( grossism * ffdx - pi ) / ( ffdx - pi )122 123 IF( 2.*beta - grossism.LE. 0.) THEN124 WRITE(6,*) ' ** Attention ! La valeur beta calculee dans la rou125 ,tine fxhyp est mauvaise ! '126 WRITE(6,*)'Modifier les valeurs de grossismx ,tau ou dzoomx ',127 , ' et relancer ! *** '128 CALL ABORT129 ENDIF130 c131 c ..... calcul de Xprimt .....132 c133 134 DO i = nmax, nmax2135 Xprimt(i) = beta + ( grossism - beta ) * fhyp(i)136 ENDDO137 c138 DO i = nmax+1, nmax2139 Xprimt( nmax2 - i ) = Xprimt( i )140 ENDDO141 c142 143 c ..... Calcul de Xf ........144 145 Xf(0) = - pi146 147 DO i = nmax +1, nmax2148 149 xmoy = 0.5 * ( xtild(i-1) + xtild( i ) )150 fa = tau* ( dzoom/2. - xmoy )151 fb = xmoy * ( pi - xmoy )152 153 IF( 200.* fb .LT. - fa ) THEN154 fxm = - 1.155 ELSEIF( 200. * fb .LT. fa ) THEN156 fxm = 1.157 ELSE158 fxm = TANH ( fa/fb )159 ENDIF160 161 192 IF ( xmoy.EQ. 0. ) fxm = 1. 162 193 IF ( xmoy.EQ. pi ) fxm = -1. … … 185 216 186 217 IF( ik.EQ.1 ) THEN 187 xuv = -0.25218 xuv = -0.25 188 219 ELSE IF ( ik.EQ.2 ) THEN 189 220 xuv = 0. 190 221 ELSE IF ( ik.EQ.3 ) THEN 191 xuv = 0.5 222 xuv = 0.50 192 223 ELSE IF ( ik.EQ.4 ) THEN 193 224 xuv = 0.25 … … 196 227 xo1 = 0. 197 228 198 DO 1500 i = 1, iim 199 200 xlon2 = - pi + (FLOAT(i) + xuv - 0.75) * depi / FLOAT(iim) 229 ii1=1 230 ii2=iim 231 IF(ik.EQ.1.and.grossism.EQ.1.) THEN 232 ii1 = 2 233 ii2 = iim+1 234 ENDIF 235 236 DO 1500 i = ii1, ii2 237 238 xlon2 = - pi + (FLOAT(i) + xuv - decalx) * depi / FLOAT(iim) 201 239 202 240 Xfi = xlon2 … … 248 286 1500 CONTINUE 249 287 288 289 IF(ik.EQ.1.and.grossism.EQ.1.) THEN 290 xvrai(1) = xvrai(iip1)-depi 291 xxprim(1) = xxprim(iip1) 292 ENDIF 250 293 DO i = 1 , iim 251 294 xlon(i) = xvrai(i) 252 295 xprimm(i) = xxprim(i) 253 296 ENDDO 254 297 255 298 DO i = 1, iim -1 256 299 IF( xvrai(i+1). LT. xvrai(i) ) THEN … … 271 314 ENDDO 272 315 273 IF(champmin .GE. - pi .AND. champmax.LE. pi) THEN316 IF(champmin .GE.-pi-0.10.and.champmax.LE.pi+0.10 ) THEN 274 317 GO TO 1600 275 318 ELSE … … 301 344 IF( ik.EQ.1 ) THEN 302 345 DO i = iim,1,-1 303 IF( xvrai(i).LE. pi )GO TO 90346 IF( xvrai(i).LE. pi ) GO TO 90 304 347 ENDDO 305 348 WRITE(6,*) ' PBS. 2 ! Xvrai plus grand que pi ! ' … … 336 379 c WRITE(6,18) 337 380 c WRITE(6,68) xvrai 338 c ccWRITE(6,*) ' XPRIM k ',ik339 c ccWRITE(6,566) xprimm340 341 DO i = 1,iim + 381 c WRITE(6,*) ' XPRIM k ',ik 382 c WRITE(6,566) xprimm 383 384 DO i = 1,iim +1 342 385 rlonm025(i) = xlon( i ) 343 386 xprimm025(i) = xprimm(i) 344 387 ENDDO 345 388 346 389 ELSE IF( ik.EQ.2 ) THEN 347 390 c WRITE(6,18) 348 391 c WRITE(6,*) ' XLON aux pts. V apres ( en deg. ) ' 349 392 c WRITE(6,68) xvrai 350 c ccWRITE(6,*) ' XPRIM k ',ik351 c ccWRITE(6,566) xprimm393 c WRITE(6,*) ' XPRIM k ',ik 394 c WRITE(6,566) xprimm 352 395 353 396 DO i = 1,iim + 1 … … 360 403 c WRITE(6,*) ' XLON aux pts. U apres ( en deg. ) ' 361 404 c WRITE(6,68) xvrai 362 c ccWRITE(6,*) ' XPRIM ik ',ik363 c ccWRITE(6,566) xprimm405 c WRITE(6,*) ' XPRIM ik ',ik 406 c WRITE(6,566) xprimm 364 407 365 408 DO i = 1,iim + 1 … … 372 415 c WRITE(6,*) ' XLON aux pts. V+0.25 apres ( en deg. ) ' 373 416 c WRITE(6,68) xvrai 374 c ccWRITE(6,*) ' XPRIM ik ',ik375 c ccWRITE(6,566) xprimm417 c WRITE(6,*) ' XPRIM ik ',ik 418 c WRITE(6,566) xprimm 376 419 377 420 DO i = 1,iim + 1 … … 403 446 24 FORMAT(2x,'Parametres xzoom,gross,tau ,dzoom pour fxhyp ',4f8.3) 404 447 68 FORMAT(1x,7f9.2) 448 566 FORMAT(1x,7f9.4) 405 449 406 450 RETURN
Note: See TracChangeset
for help on using the changeset viewer.
