/[lmdze]/trunk/Sources/dyn3d/invert_zoom_x.f
ViewVC logotype

Annotation of /trunk/Sources/dyn3d/invert_zoom_x.f

Parent Directory Parent Directory | Revision Log Revision Log

Revision 144 - (hide annotations)
Wed Jun 10 16:46:46 2015 UTC (8 years, 11 months ago) by guez
File size: 2208 byte(s) 
In procedure fxhyp, the convoluted computation of tanh(fa/fb) occurred
three times. Extracted it into a function. Also, the computation of
xmoy and fxm was repeated. So stored the values in arrays instead.

In procedure fxhyp, in the computation of fhyp, there were tests
xtild(i) == 0. and xtild(i) == pi_d. No use to do these tests at each
iteration. We now they are true for i == nmax and i == 2 * nmax,
respectively, and we know they are false for other values of
"i". Similarly, in the computations of ffdx and xxpr, there were the
tests xmoy == 0. and xmoy == pi_d, these could not be true.

Moved files from bibio to dyn3d, following LMDZ.
1 guez 131 module invert_zoom_x_m
2 guez 121
3     implicit none
4    
5     INTEGER, PARAMETER:: nmax = 30000
6    
7     contains
8    
9 guez 131 subroutine invert_zoom_x(xf, xtild, Xprimt, xlon, xprimm, xuv)
10 guez 121
11     use coefpoly_m, only: coefpoly
12     USE dimens_m, ONLY: iim
13 guez 139 use dynetat0_m, only: clon
14 guez 124 use nr_util, only: pi_d, twopi_d
15 guez 121
16     DOUBLE PRECISION, intent(in):: Xf(0:), xtild(0:), Xprimt(0:) ! (0:2 * nmax)
17 guez 124 real, intent(out):: xlon(:), xprimm(:) ! (iim)
18 guez 121
19     DOUBLE PRECISION, intent(in):: xuv
20     ! 0. si calcul aux points scalaires
21     ! 0.5 si calcul aux points U
22    
23     ! Local:
24 guez 126 DOUBLE PRECISION xo1, Xfi, a0, a1, a2, a3, Xf1, Xprimin
25     integer i, it, iter
26     DOUBLE PRECISION, parameter:: my_eps = 1d-6
27 guez 121
28 guez 126 DOUBLE PRECISION xxprim(iim), xvrai(iim)
29     ! intermediary variables because xlon and xprimm are simple precision
30    
31 guez 121 !------------------------------------------------------------------
32    
33 guez 126 DO i = 1, iim
34     Xfi = - pi_d + (i + xuv - 0.75d0) * twopi_d / iim
35 guez 121
36     it = 2 * nmax
37     do while (xfi < xf(it) .and. it >= 1)
38     it = it - 1
39     end do
40    
41 guez 126 ! Calcul de Xf(xvrai(i))
42 guez 121
43 guez 126 xvrai(i) = xtild(it)
44 guez 144 IF (it == 2 * nmax) it = 2 * nmax -1
45 guez 121 CALL coefpoly(Xf(it), Xf(it + 1), Xprimt(it), Xprimt(it + 1), &
46     xtild(it), xtild(it + 1), a0, a1, a2, a3)
47     Xf1 = Xf(it)
48 guez 126 Xprimin = a1 + xvrai(i) * (2d0 * a2 + xvrai(i) * 3d0 * a3)
49     xo1 = xvrai(i)
50 guez 121 iter = 1
51    
52     do
53 guez 126 xvrai(i) = xvrai(i) - (Xf1 - Xfi) / Xprimin
54     IF (ABS(xvrai(i) - xo1) <= my_eps .or. iter == 300) exit
55     xo1 = xvrai(i)
56     Xf1 = a0 + xvrai(i) * (a1 + xvrai(i) * (a2 + xvrai(i) * a3))
57     Xprimin = a1 + xvrai(i) * (2d0 * a2 + xvrai(i) * 3d0 * a3)
58 guez 121 end DO
59    
60 guez 126 if (ABS(xvrai(i) - xo1) > my_eps) then
61 guez 121 ! iter == 300
62     print *, 'Pas de solution.'
63     print *, i, xfi
64     STOP 1
65     end if
66    
67 guez 126 xxprim(i) = twopi_d / (iim * Xprimin)
68 guez 121 end DO
69    
70     DO i = 1, iim -1
71     IF (xvrai(i + 1) < xvrai(i)) THEN
72 guez 124 print *, 'xvrai(', i + 1, ') < xvrai(', i, ')'
73 guez 121 STOP 1
74     END IF
75     END DO
76    
77 guez 127 xlon = xvrai + clon
78 guez 126 xprimm = xxprim
79 guez 121
80 guez 131 end subroutine invert_zoom_x
81 guez 121
82 guez 131 end module invert_zoom_x_m
  ViewVC Help
Powered by ViewVC 1.1.21