source: codes/icosagcm/trunk/src/domain.f90 @ 21

Last change on this file since 21 was 21, checked in by ymipsl, 12 years ago

correction for compiling with gfortran (line too long)
improvement for splitting domain
Call twice transfert request for field u is no longer necessary

YM
File size: 8.5 KB
Line 
1MODULE domain_mod
2  USE domain_param
3
4  TYPE t_domain
5    INTEGER :: face
6    INTEGER :: iim
7    INTEGER :: jjm
8    INTEGER :: ii_begin
9    INTEGER :: jj_begin
10    INTEGER :: ii_end
11    INTEGER :: jj_end
12    INTEGER :: ii_nb
13    INTEGER :: jj_nb
14    INTEGER :: ii_begin_glo
15    INTEGER :: jj_begin_glo
16    INTEGER :: ii_end_glo
17    INTEGER :: jj_end_glo
18    INTEGER,POINTER  :: assign_domain(:,:)
19    INTEGER,POINTER  :: assign_i(:,:)
20    INTEGER,POINTER  :: assign_j(:,:)
21    INTEGER,POINTER  :: edge_assign_domain(:,:,:)
22    INTEGER,POINTER  :: edge_assign_i(:,:,:)
23    INTEGER,POINTER  :: edge_assign_j(:,:,:)
24    INTEGER,POINTER  :: edge_assign_pos(:,:,:)
25    REAL,POINTER     :: xyz(:,:,:)
26    REAL,POINTER     :: neighbour(:,:,:,:)
27    INTEGER,POINTER  :: delta(:,:)
28    REAL,POINTER     :: vertex(:,:,:,:)
29    LOGICAL,POINTER  :: own(:,:)
30    INTEGER,POINTER  :: ne(:,:,:)
31   
32    INTEGER :: t_right
33    INTEGER :: t_rup
34    INTEGER :: t_lup
35    INTEGER :: t_left
36    INTEGER :: t_ldown
37    INTEGER :: t_rdown
38
39    INTEGER :: u_right
40    INTEGER :: u_rup
41    INTEGER :: u_lup
42    INTEGER :: u_left
43    INTEGER :: u_ldown
44    INTEGER :: u_rdown
45
46    INTEGER :: z_rup
47    INTEGER :: z_up
48    INTEGER :: z_lup
49    INTEGER :: z_ldown
50    INTEGER :: z_down
51    INTEGER :: z_rdown
52   
53    INTEGER :: t_pos(6)
54    INTEGER :: u_pos(6)
55    INTEGER :: z_pos(6)
56     
57  END TYPE t_domain
58
59  TYPE(t_domain),SAVE,ALLOCATABLE,TARGET :: domain(:)
60  INTEGER :: ndomain
61
62
63CONTAINS
64
65  SUBROUTINE create_domain
66  USE grid_param
67  IMPLICIT NONE
68  INTEGER :: ind,nf,ni,nj,i,j
69  INTEGER :: quotient, rest
70  TYPE(t_domain),POINTER :: d
71 
72    ndomain=nsplit_i*nsplit_j*nb_face
73    ALLOCATE(domain(ndomain))
74 
75    ind=0
76    DO nf=1,nb_face
77      DO nj=1,nsplit_j
78        DO ni=1,nsplit_i
79          ind=ind+1
80          d=>domain(ind)
81          quotient=(iim_glo/nsplit_i)
82          rest=MOD(iim_glo,nsplit_i)
83          IF (ni-1 < rest) THEN
84            d%ii_nb=quotient+1
85            d%ii_begin_glo=(quotient+1)*(ni-1)+1
86          ELSE
87            d%ii_nb=quotient
88            d%ii_begin_glo=(quotient+1)*rest+(ni-1-rest) * quotient+1
89          ENDIF
90          d%ii_end_glo=d%ii_begin_glo+d%ii_nb-1
91 
92          IF (ni/=nsplit_i) THEN
93            d%ii_nb=d%ii_nb+1
94            d%ii_end_glo=d%ii_end_glo+1
95          ENDIF
96         
97       
98          quotient=(jjm_glo/nsplit_j)
99          rest=MOD(jjm_glo,nsplit_j)
100          IF (nj-1 < rest) THEN
101            d%jj_nb=quotient+1
102            d%jj_begin_glo=(quotient+1)*(nj-1)+1
103          ELSE
104            d%jj_nb=quotient
105            d%jj_begin_glo=(quotient+1)*rest+(nj-1-rest) * quotient+1
106          ENDIF
107          d%jj_end_glo=d%jj_begin_glo+d%jj_nb-1
108
109          IF (nj/=nsplit_j) THEN
110            d%jj_nb=d%jj_nb+1
111            d%jj_end_glo=d%jj_end_glo+1
112          ENDIF
113
114
115          d%iim=d%ii_nb+2*halo
116          d%jjm=d%jj_nb+2*halo
117          d%ii_begin=halo+1
118          d%jj_begin=halo+1
119          d%ii_end=d%ii_begin+d%ii_nb-1
120          d%jj_end=d%jj_begin+d%jj_nb-1
121          d%face=nf       
122          ALLOCATE(d%assign_domain(d%iim,d%jjm))
123          ALLOCATE(d%assign_i(d%iim,d%jjm))
124          ALLOCATE(d%assign_j(d%iim,d%jjm))
125          ALLOCATE(d%edge_assign_domain(0:5,d%iim,d%jjm))
126          ALLOCATE(d%edge_assign_i(0:5,d%iim,d%jjm))
127          ALLOCATE(d%edge_assign_j(0:5,d%iim,d%jjm))
128          ALLOCATE(d%edge_assign_pos(0:5,d%iim,d%jjm))
129          ALLOCATE(d%delta(d%iim,d%jjm))
130          ALLOCATE(d%xyz(3,d%iim,d%jjm))
131          ALLOCATE(d%vertex(3,0:5,d%iim,d%jjm))
132          ALLOCATE(d%neighbour(3,0:5,d%iim,d%jjm))
133          ALLOCATE(d%own(d%iim,d%jjm))
134          ALLOCATE(d%ne(0:5,d%iim,d%jjm))
135        END DO
136      END DO
137    END DO
138   
139  END SUBROUTINE create_domain
140 
141 
142  SUBROUTINE assign_cell
143  USE metric
144  IMPLICIT NONE
145    INTEGER :: ind_d,ind,ind2,e
146    INTEGER :: nf
147    INTEGER :: i,j,k,ii,jj
148    TYPE(t_domain),POINTER :: d
149    INTEGER :: delta
150     
151   
152    DO ind_d=1,ndomain
153      d=>domain(ind_d)
154      nf=d%face
155      DO j=d%jj_begin,d%jj_end
156        DO i=d%ii_begin,d%ii_end
157          ii=d%ii_begin_glo-d%ii_begin+i
158          jj=d%jj_begin_glo-d%jj_begin+j
159          ind=vertex_glo(ii,jj,nf)%ind
160          delta=vertex_glo(ii,jj,nf)%delta
161          IF (cell_glo(ind)%assign_face==nf) THEN
162            cell_glo(ind)%assign_domain=ind_d
163            cell_glo(ind)%assign_i=i
164            cell_glo(ind)%assign_j=j
165          ENDIF
166          IF ( i+1 <= d%ii_end ) CALL assign_edge(ind_d,ind,i,j,delta,0)
167          IF ( j+1 <= d%jj_end ) CALL assign_edge(ind_d,ind,i,j,delta,1)
168          IF ( i-1 >= d%ii_begin .AND. j+1<=d%jj_end ) CALL assign_edge(ind_d,ind,i,j,delta,2)
169          IF ( i-1 >= d%ii_begin ) CALL assign_edge(ind_d,ind,i,j,delta,3)
170          IF ( j-1 >= d%jj_begin ) CALL assign_edge(ind_d,ind,i,j,delta,4)
171          IF ( i+1 <= d%ii_end .AND. j-1 >=d%jj_begin ) CALL assign_edge(ind_d,ind,i,j,delta,5)
172        ENDDO
173      ENDDO
174    ENDDO
175   
176   
177    DO ind_d=1,ndomain
178      d=>domain(ind_d)
179      nf=d%face
180      DO j=d%jj_begin-1,d%jj_end+1
181        DO i=d%ii_begin-1,d%ii_end+1
182          ii=d%ii_begin_glo-d%ii_begin+i
183          jj=d%jj_begin_glo-d%jj_begin+j
184          ind=vertex_glo(ii,jj,nf)%ind
185          d%assign_domain(i,j)=cell_glo(ind)%assign_domain
186          d%assign_i(i,j)=cell_glo(ind)%assign_i
187          d%assign_j(i,j)=cell_glo(ind)%assign_j
188          delta=vertex_glo(ii,jj,nf)%delta
189          d%delta(i,j)=vertex_glo(ii,jj,nf)%delta
190          DO k=0,5
191            ind2=vertex_glo(ii,jj,nf)%neighbour(k)
192            d%neighbour(:,k,i,j)=cell_glo(ind2)%xyz(:)
193            d%ne(k,i,j)=vertex_glo(ii,jj,nf)%ne(k)
194            e=cell_glo(ind)%edge(MOD(k+delta+6,6))
195            d%edge_assign_domain(k,i,j)=edge_glo(e)%assign_domain
196            d%edge_assign_i(k,i,j)=edge_glo(e)%assign_i
197            d%edge_assign_j(k,i,j)=edge_glo(e)%assign_j
198            d%edge_assign_pos(k,i,j)=edge_glo(e)%assign_pos
199          ENDDO
200          d%xyz(:,i,j)=cell_glo(ind)%xyz(:)
201          IF (d%assign_domain(i,j)==ind_d) THEN
202           d%own(i,j)=.TRUE.
203          ELSE
204           d%own(i,j)=.FALSE.
205          ENDIF
206        ENDDO
207      ENDDO
208    ENDDO
209
210  CONTAINS
211
212    SUBROUTINE assign_edge(ind_d,ind,i,j,delta,k)
213      INTEGER :: ind_d,ind,i,j,delta,k
214      INTEGER :: e
215      e=cell_glo(ind)%edge(MOD(k+delta+6,6))
216      edge_glo(e)%assign_domain=ind_d
217      edge_glo(e)%assign_i=i
218      edge_glo(e)%assign_j=j
219      edge_glo(e)%assign_pos=k
220!      PRINT*,"Assign_edge",ind_d,ind,i,j,k,e
221     END  SUBROUTINE assign_edge
222         
223  END SUBROUTINE assign_cell
224
225  SUBROUTINE compute_boundary
226  USE spherical_geom_mod
227  IMPLICIT NONE
228    INTEGER :: ind_d
229    INTEGER :: i,j,k
230    TYPE(t_domain),POINTER :: d 
231    REAL(rstd) :: ng1(3),ng2(3) 
232
233    DO ind_d=1,ndomain
234      d=>domain(ind_d)
235      DO j=d%jj_begin-1,d%jj_end+1
236        DO i=d%ii_begin-1,d%ii_end+1
237          DO k=0,5
238            ng1=d%neighbour(:,MOD(k,6),i,j)
239            ng2=d%neighbour(:,MOD(k+1,6),i,j)
240            IF (sqrt(sum((ng1-ng2)**2))<1e-16) ng2=d%neighbour(:,MOD(k+2,6),i,j)
241            CALL circumcenter(d%xyz(:,i,j),ng1,ng2,d%vertex(:,k,i,j))
242          ENDDO
243        ENDDO
244      ENDDO
245    ENDDO       
246  END SUBROUTINE compute_boundary
247
248  SUBROUTINE set_neighbour_indice
249  USE metric
250  IMPLICIT NONE
251    INTEGER :: ind_d
252    TYPE(t_domain),POINTER :: d 
253   
254    DO ind_d=1,ndomain
255      d=>domain(ind_d)
256      d%t_right=1
257      d%t_left=-1
258      d%t_rup=d%iim
259      d%t_lup=d%iim-1
260      d%t_ldown=-d%iim
261      d%t_rdown=-d%iim+1
262     
263      d%u_right=0
264      d%u_lup=d%iim*d%jjm
265      d%u_ldown=2*d%iim*d%jjm
266     
267      d%u_rup=d%t_rup+d%u_ldown
268      d%u_left=d%t_left+d%u_right
269      d%u_rdown=d%t_rdown+d%u_lup
270     
271      d%z_up=0
272      d%z_down=d%iim*d%jjm
273      d%z_rup=d%t_rup+d%z_down
274      d%z_lup=d%t_lup+d%z_down
275      d%z_ldown=d%t_ldown+d%z_up
276      d%z_rdown=d%t_rdown+d%z_up
277     
278      d%t_pos(right)=d%t_right
279      d%t_pos(rup)=d%t_rup
280      d%t_pos(lup)=d%t_lup
281      d%t_pos(left)=d%t_left
282      d%t_pos(ldown)=d%t_ldown
283      d%t_pos(rdown)=d%t_rdown
284
285      d%u_pos(right)=d%u_right
286      d%u_pos(rup)=d%u_rup
287      d%u_pos(lup)=d%u_lup
288      d%u_pos(left)=d%u_left
289      d%u_pos(ldown)=d%u_ldown
290      d%u_pos(rdown)=d%u_rdown
291     
292      d%z_pos(vrup)=d%z_rup
293      d%z_pos(vup)=d%z_up
294      d%z_pos(vlup)=d%z_lup
295      d%z_pos(vldown)=d%z_ldown
296      d%z_pos(vdown)=d%z_down
297      d%z_pos(vrdown)=d%z_rdown
298     
299    ENDDO 
300 
301  END SUBROUTINE set_neighbour_indice
302     
303     
304         
305  SUBROUTINE compute_domain
306  IMPLICIT NONE
307 
308    CALL create_domain
309    CALL assign_cell
310    CALL compute_boundary
311    CALL set_neighbour_indice
312     
313  END SUBROUTINE compute_domain
314         
315END MODULE domain_mod 
316 
Note: See TracBrowser for help on using the repository browser.