My Project
 All Classes Files Functions Variables Macros
mod_filtre_fft_loc.F90
Go to the documentation of this file.
1 MODULE mod_filtre_fft_loc
2 
3  LOGICAL,SAVE :: use_filtre_fft
4  REAL,SAVE,ALLOCATABLE :: Filtre_u(:,:)
5  REAL,SAVE,ALLOCATABLE :: Filtre_v(:,:)
6  REAL,SAVE,ALLOCATABLE :: Filtre_inv(:,:)
7 
8 CONTAINS
9 
10  SUBROUTINE init_filtre_fft(coeffu,modfrstu,jfiltnu,jfiltsu,coeffv,modfrstv,jfiltnv,jfiltsv)
11  USE mod_fft
12  IMPLICIT NONE
13  include 'dimensions.h'
14  REAL, INTENT(IN) :: coeffu(iim,jjm)
15  INTEGER,INTENT(IN) :: modfrstu(jjm)
16  INTEGER,INTENT(IN) :: jfiltnu
17  INTEGER,INTENT(IN) :: jfiltsu
18  REAL, INTENT(IN) :: coeffv(iim,jjm)
19  INTEGER,INTENT(IN) :: modfrstv(jjm)
20  INTEGER,INTENT(IN) :: jfiltnv
21  INTEGER,INTENT(IN) :: jfiltsv
22 
23  INTEGER :: index_vp(iim)
24  INTEGER :: i,j
25 
26  index_vp(1)=1
27  DO i=1,iim/2
28  index_vp(i+1)=i*2
29  ENDDO
30 
31  DO i=1,iim/2-1
32  index_vp(iim/2+i+1)=iim-2*i+1
33  ENDDO
34 
35  ALLOCATE(filtre_u(iim,jjm))
36  ALLOCATE(filtre_v(iim,jjm))
37  ALLOCATE(filtre_inv(iim,jjm))
38 
39 
40  DO j=2,jfiltnu
41  DO i=1,iim
42  IF (index_vp(i) < modfrstu(j)) THEN
43  filtre_u(i,j)=0
44  ELSE
45  filtre_u(i,j)=coeffu(index_vp(i),j)
46  ENDIF
47  ENDDO
48  ENDDO
49 
50  DO j=jfiltsu,jjm
51  DO i=1,iim
52  IF (index_vp(i) < modfrstu(j)) THEN
53  filtre_u(i,j)=0
54  ELSE
55  filtre_u(i,j)=coeffu(index_vp(i),j)
56  ENDIF
57  ENDDO
58  ENDDO
59 
60  DO j=1,jfiltnv
61  DO i=1,iim
62  IF (index_vp(i) < modfrstv(j)) THEN
63  filtre_v(i,j)=0
64  ELSE
65  filtre_v(i,j)=coeffv(index_vp(i),j)
66  ENDIF
67  ENDDO
68  ENDDO
69 
70  DO j=jfiltsv,jjm
71  DO i=1,iim
72  IF (index_vp(i) < modfrstv(j)) THEN
73  filtre_v(i,j)=0
74  ELSE
75  filtre_v(i,j)=coeffv(index_vp(i),j)
76  ENDIF
77  ENDDO
78  ENDDO
79 
80  DO j=2,jfiltnu
81  DO i=1,iim
82  IF (index_vp(i) < modfrstu(j)) THEN
83  filtre_inv(i,j)=0
84  ELSE
85  filtre_inv(i,j)=coeffu(index_vp(i),j)/(1.+coeffu(index_vp(i),j))
86  ENDIF
87  ENDDO
88  ENDDO
89 
90  DO j=jfiltsu,jjm
91  DO i=1,iim
92  IF (index_vp(i) < modfrstu(j)) THEN
93  filtre_inv(i,j)=0
94  ELSE
95  filtre_inv(i,j)=coeffu(index_vp(i),j)/(1.+coeffu(index_vp(i),j))
96  ENDIF
97  ENDDO
98  ENDDO
99 
100 
101 ! CALL Init_FFT(iim,(jjm+1)*(llm+1))
102 
103 
104  END SUBROUTINE init_filtre_fft
105 
106  SUBROUTINE filtre_u_fft(vect_inout,jjb,jje,jj_begin,jj_end,nbniv)
107  USE mod_fft
108 #ifdef CPP_PARA
109  USE parallel,ONLY : OMP_CHUNK
110 #endif
111  IMPLICIT NONE
112  include 'dimensions.h'
113  INTEGER,INTENT(IN) :: jjb
114  INTEGER,INTENT(IN) :: jje
115  INTEGER,INTENT(IN) :: jj_begin
116  INTEGER,INTENT(IN) :: jj_end
117  INTEGER,INTENT(IN) :: nbniv
118  REAL,INTENT(INOUT) :: vect_inout(iim+1,jjb:jje,nbniv)
119 
120  REAL :: vect(iim+inc,jj_end-jj_begin+1,nbniv)
121  COMPLEX :: tf_vect(iim/2+1,jj_end-jj_begin+1,nbniv)
122  INTEGER :: nb_vect
123  INTEGER :: i,j,l
124  INTEGER :: ll_nb
125 ! REAL :: vect_tmp(iim+inc,jj_end-jj_begin+1,nbniv)
126 
127  ll_nb=0
128 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
129  DO l=1,nbniv
130  ll_nb=ll_nb+1
131  DO j=1,jj_end-jj_begin+1
132  DO i=1,iim+1
133  vect(i,j,ll_nb)=vect_inout(i,j+jj_begin-1,l)
134  ENDDO
135  ENDDO
136  ENDDO
137 !$OMP END DO NOWAIT
138 
139  nb_vect=(jj_end-jj_begin+1)*ll_nb
140 
141 ! vect_tmp=vect
142 
143  CALL fft_forward(vect,tf_vect,nb_vect)
144 
145 ! CALL FFT_forward(vect,TF_vect_test,nb_vect)
146 ! PRINT *,"XXXXXXXXXXXXX Filtre_u_FFT xxxxxxxxxxxx"
147 ! DO j=1,jj_end-jj_begin+1
148 ! DO i=1,iim/2+1
149 ! PRINT *,"====",i,j,"----->",TF_vect_test(i,j,1)
150 ! ENDDO
151 ! ENDDO
152 
153  DO l=1,ll_nb
154  DO j=1,jj_end-jj_begin+1
155  DO i=1,iim/2+1
156  tf_vect(i,j,l)=tf_vect(i,j,l)*filtre_u(i,jj_begin+j-1)
157  ENDDO
158  ENDDO
159  ENDDO
160 
161  CALL fft_backward(tf_vect,vect,nb_vect)
162 ! CALL FFT_backward(TF_vect_test,vect_test,nb_vect)
163 
164 ! PRINT *,"XXXXXXXXXXXXX Filtre_u_FFT xxxxxxxxxxxx"
165 ! DO j=1,jj_end-jj_begin+1
166 ! DO i=1,iim
167 ! PRINT *,"====",i,j,"----->",vect_test(i,j,1)
168 ! ENDDO
169 ! ENDDO
170 
171  ll_nb=0
172 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
173  DO l=1,nbniv
174  ll_nb=ll_nb+1
175  DO j=1,jj_end-jj_begin+1
176  DO i=1,iim+1
177  vect_inout(i,j+jj_begin-1,l)=vect(i,j,ll_nb)
178  ENDDO
179  ENDDO
180  ENDDO
181 !$OMP END DO NOWAIT
182 
183  END SUBROUTINE filtre_u_fft
184 
185 
186  SUBROUTINE filtre_v_fft(vect_inout,jjb,jje,jj_begin,jj_end,nbniv)
187  USE mod_fft
188 #ifdef CPP_PARA
189  USE parallel,ONLY : OMP_CHUNK
190 #endif
191  IMPLICIT NONE
192  include 'dimensions.h'
193  INTEGER,INTENT(IN) :: jjb
194  INTEGER,INTENT(IN) :: jje
195  INTEGER,INTENT(IN) :: jj_begin
196  INTEGER,INTENT(IN) :: jj_end
197  INTEGER,INTENT(IN) :: nbniv
198  REAL,INTENT(INOUT) :: vect_inout(iim+1,jjb:jje,nbniv)
199 
200  REAL :: vect(iim+inc,jj_end-jj_begin+1,nbniv)
201  COMPLEX :: tf_vect(iim/2+1,jj_end-jj_begin+1,nbniv)
202  INTEGER :: nb_vect
203  INTEGER :: i,j,l
204  INTEGER :: ll_nb
205 
206  ll_nb=0
207 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
208  DO l=1,nbniv
209  ll_nb=ll_nb+1
210  DO j=1,jj_end-jj_begin+1
211  DO i=1,iim+1
212  vect(i,j,ll_nb)=vect_inout(i,j+jj_begin-1,l)
213  ENDDO
214  ENDDO
215  ENDDO
216 !$OMP END DO NOWAIT
217 
218 
219  nb_vect=(jj_end-jj_begin+1)*ll_nb
220 
221  CALL fft_forward(vect,tf_vect,nb_vect)
222 
223  DO l=1,ll_nb
224  DO j=1,jj_end-jj_begin+1
225  DO i=1,iim/2+1
226  tf_vect(i,j,l)=tf_vect(i,j,l)*filtre_v(i,jj_begin+j-1)
227  ENDDO
228  ENDDO
229  ENDDO
230 
231  CALL fft_backward(tf_vect,vect,nb_vect)
232 
233 
234  ll_nb=0
235 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
236  DO l=1,nbniv
237  ll_nb=ll_nb+1
238  DO j=1,jj_end-jj_begin+1
239  DO i=1,iim+1
240  vect_inout(i,j+jj_begin-1,l)=vect(i,j,ll_nb)
241  ENDDO
242  ENDDO
243  ENDDO
244 !$OMP END DO NOWAIT
245 
246  END SUBROUTINE filtre_v_fft
247 
248 
249  SUBROUTINE filtre_inv_fft(vect_inout,jjb,jje,jj_begin,jj_end,nbniv)
250  USE mod_fft
251 #ifdef CPP_PARA
252  USE parallel,ONLY : OMP_CHUNK
253 #endif
254  IMPLICIT NONE
255  include 'dimensions.h'
256  INTEGER,INTENT(IN) :: jjb
257  INTEGER,INTENT(IN) :: jje
258  INTEGER,INTENT(IN) :: jj_begin
259  INTEGER,INTENT(IN) :: jj_end
260  INTEGER,INTENT(IN) :: nbniv
261  REAL,INTENT(INOUT) :: vect_inout(iim+1,jjb:jje,nbniv)
262 
263  REAL :: vect(iim+inc,jj_end-jj_begin+1,nbniv)
264  COMPLEX :: tf_vect(iim/2+1,jj_end-jj_begin+1,nbniv)
265  INTEGER :: nb_vect
266  INTEGER :: i,j,l
267  INTEGER :: ll_nb
268 
269  ll_nb=0
270 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
271  DO l=1,nbniv
272  ll_nb=ll_nb+1
273  DO j=1,jj_end-jj_begin+1
274  DO i=1,iim+1
275  vect(i,j,ll_nb)=vect_inout(i,j+jj_begin-1,l)
276  ENDDO
277  ENDDO
278  ENDDO
279 !$OMP END DO NOWAIT
280 
281  nb_vect=(jj_end-jj_begin+1)*ll_nb
282 
283  CALL fft_forward(vect,tf_vect,nb_vect)
284 
285  DO l=1,ll_nb
286  DO j=1,jj_end-jj_begin+1
287  DO i=1,iim/2+1
288  tf_vect(i,j,l)=tf_vect(i,j,l)*filtre_inv(i,jj_begin+j-1)
289  ENDDO
290  ENDDO
291  ENDDO
292 
293  CALL fft_backward(tf_vect,vect,nb_vect)
294 
295  ll_nb=0
296 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
297  DO l=1,nbniv
298  ll_nb=ll_nb+1
299  DO j=1,jj_end-jj_begin+1
300  DO i=1,iim+1
301  vect_inout(i,j+jj_begin-1,l)=vect(i,j,ll_nb)
302  ENDDO
303  ENDDO
304  ENDDO
305 !$OMP END DO NOWAIT
306 
307  END SUBROUTINE filtre_inv_fft
308 
309 
310 ! SUBROUTINE get_ll_index(nbniv,ll_index,ll_nb)
311 ! IMPLICIT NONE
312 ! INTEGER,INTENT(IN) :: nbniv
313 ! INTEGER,INTENT(OUT) :: ll_index(nbniv)
314 ! INTEGER,INTENT(OUT) :: ll_nb
315 !
316 ! INTEGER :: l,ll_begin, ll_end
317 ! INTEGER :: omp_rank,omp_size
318 ! INTEGER :: OMP_GET_NUM_THREADS
319 ! INTEGER :: omp_chunk
320 ! EXTERNAL OMP_GET_NUM_THREADS
321 ! INTEGER :: OMP_GET_THREAD_NUM
322 ! EXTERNAL OMP_GET_THREAD_NUM
323 !
324 !
325 ! omp_size=OMP_GET_NUM_THREADS()
326 ! omp_rank=OMP_GET_THREAD_NUM()
327 ! omp_chunk=nbniv/omp_size+min(1,MOD(nbniv,omp_size))
328 !
329 ! ll_begin=omp_rank*OMP_CHUNK+1
330 ! ll_nb=0
331 ! DO WHILE (ll_begin<=nbniv)
332 ! ll_end=min(ll_begin+OMP_CHUNK-1,nbniv)
333 ! DO l=ll_begin,ll_end
334 ! ll_nb=ll_nb+1
335 ! ll_index(ll_nb)=l
336 ! ENDDO
337 ! ll_begin=ll_begin+omp_size*OMP_CHUNK
338 ! ENDDO
339 !
340 ! END SUBROUTINE get_ll_index
341 
342 END MODULE mod_filtre_fft_loc
343