/[lmdze]/trunk/Sources/dyn3d/PPM3d/filns.f
ViewVC logotype

Contents of /trunk/Sources/dyn3d/PPM3d/filns.f

Parent Directory Parent Directory | Revision Log Revision Log

Revision 166 - (show annotations)
Wed Jul 29 14:32:55 2015 UTC (8 years, 10 months ago) by guez
File size: 1759 byte(s) 
Split ppm3d.f into files containing a single procedure.

Factorized computations of filtering matrices into a procedure
inifilr_hemisph. Had then to change the matrices from allocatable to
pointer and from customized lower bound to lower bound 1. The change
in lower bounds does not matter because the matrices are only used as
a whole as actual arguments.

Also, in infilr, instead of finding jfilt[ns][uv] from approximately
jjm /2, start at index j1 that corresponds to the equator. This is not
the same if there is a zoom in latitude.

Also, the test (rlamda(modfrst[ns][uv](j)) * cos(rlat[uv](j)) < 1) in
the loops on filtered latitudes is not useful now that we start from
j1: it is necessarily true. See notes.

Just encapsulated lwvn into a module and removed unused argument ktraer.
1 SUBROUTINE filns(q, imr, jnp, j1, j2, cosp, acosp, ipy, tiny)
2 DIMENSION q(imr, *), cosp(*), acosp(*)
3 ! logical first
4 ! data first /.true./
5 ! save cap1
6
7 ! if(first) then
8 dp = 4.*atan(1.)/float(jnp-1)
9 cap1 = imr*(1.-cos((j1-1.5)*dp))/dp
10 ! first = .false.
11 ! endif
12
13 ipy = 0
14 DO j = j1 + 1, j2 - 1
15 DO i = 1, imr
16 IF (q(i,j)<0.) THEN
17 ipy = 1
18 dq = -q(i, j)*cosp(j)
19 ! North
20 dn = q(i, j+1)*cosp(j+1)
21 d0 = max(0., dn)
22 d1 = min(dq, d0)
23 q(i, j+1) = (dn-d1)*acosp(j+1)
24 dq = dq - d1
25 ! South
26 ds = q(i, j-1)*cosp(j-1)
27 d0 = max(0., ds)
28 d2 = min(dq, d0)
29 q(i, j-1) = (ds-d2)*acosp(j-1)
30 q(i, j) = (d2-dq)*acosp(j) + tiny
31 END IF
32 END DO
33 END DO
34
35 DO i = 1, imr
36 IF (q(i,j1)<0.) THEN
37 ipy = 1
38 dq = -q(i, j1)*cosp(j1)
39 ! North
40 dn = q(i, j1+1)*cosp(j1+1)
41 d0 = max(0., dn)
42 d1 = min(dq, d0)
43 q(i, j1+1) = (dn-d1)*acosp(j1+1)
44 q(i, j1) = (d1-dq)*acosp(j1) + tiny
45 END IF
46 END DO
47
48 j = j2
49 DO i = 1, imr
50 IF (q(i,j)<0.) THEN
51 ipy = 1
52 dq = -q(i, j)*cosp(j)
53 ! South
54 ds = q(i, j-1)*cosp(j-1)
55 d0 = max(0., ds)
56 d2 = min(dq, d0)
57 q(i, j-1) = (ds-d2)*acosp(j-1)
58 q(i, j) = (d2-dq)*acosp(j) + tiny
59 END IF
60 END DO
61
62 ! Check Poles.
63 IF (q(1,1)<0.) THEN
64 dq = q(1, 1)*cap1/float(imr)*acosp(j1)
65 DO i = 1, imr
66 q(i, 1) = 0.
67 q(i, j1) = q(i, j1) + dq
68 IF (q(i,j1)<0.) ipy = 1
69 END DO
70 END IF
71
72 IF (q(1,jnp)<0.) THEN
73 dq = q(1, jnp)*cap1/float(imr)*acosp(j2)
74 DO i = 1, imr
75 q(i, jnp) = 0.
76 q(i, j2) = q(i, j2) + dq
77 IF (q(i,j2)<0.) ipy = 1
78 END DO
79 END IF
80
81 RETURN
82 END SUBROUTINE filns
83
  ViewVC Help
Powered by ViewVC 1.1.21