Context Navigation

source: codes/icosagcm/trunk/src/metric.f90 @ 146

Last change on this file since 146 was 146, checked in by ymipsl, 11 years ago

Set constant sign for wind way :
ne(ij,right)==ne_right=1
ne(ij,rup)==ne_rup=-1
ne(ij,lup)==ne_lup=1
ne(ij,left)==ne_left=-1
ne(ij,ldown)==ne_ldown=1
ne(ij,rdown)==ne_rdown=-1

Modified transfert function to be compliant for this convention.

YM

File size: 24.6 KB

Line
1	MODULE metric
2	USE genmod
3	USE grid_param
4	USE domain_param
5
6	TYPE t_cell_glo
7	INTEGER :: neighbour(0:5)
8	INTEGER :: edge(0:5)
9	INTEGER :: assign_face
10	INTEGER :: assign_i
11	INTEGER :: assign_j
12	INTEGER :: assign_domain
13	REAL(rstd) :: xyz(3)
14	END TYPE t_cell_glo
15
16	TYPE t_vertex_glo
17	REAL(rstd) :: xyz(3)
18	INTEGER :: neighbour(0:5)
19	INTEGER :: ne(0:5)
20	INTEGER :: ind
21	INTEGER :: delta
22	END TYPE t_vertex_glo
23
24	TYPE t_edge_glo
25	INTEGER :: e1
26	INTEGER :: e2
27	INTEGER :: assign_domain
28	INTEGER :: assign_i
29	INTEGER :: assign_j
30	INTEGER :: assign_pos
31	INTEGER :: assign_delta
32	END TYPE t_edge_glo
33
34
35	TYPE(t_vertex_glo),ALLOCATABLE,SAVE :: vertex_glo(:,:,:)
36	TYPE(t_cell_glo),ALLOCATABLE,SAVE :: cell_glo(:)
37	TYPE(t_edge_glo),ALLOCATABLE,SAVE :: edge_glo(:)
38	INTEGER :: ncell_glo
39
40	INTEGER,ALLOCATABLE,SAVE :: Tab_index(:,:,:)
41	INTEGER,PARAMETER :: right=1
42	INTEGER,PARAMETER :: rup=2
43	INTEGER,PARAMETER :: lup=3
44	INTEGER,PARAMETER :: left=4
45	INTEGER,PARAMETER :: ldown=5
46	INTEGER,PARAMETER :: rdown=6
47
48	INTEGER,PARAMETER :: vrup=1
49	INTEGER,PARAMETER :: vup=2
50	INTEGER,PARAMETER :: vlup=3
51	INTEGER,PARAMETER :: vldown=4
52	INTEGER,PARAMETER :: vdown=5
53	INTEGER,PARAMETER :: vrdown=6
54
55
56	CONTAINS
57
58	SUBROUTINE allocate_metric
59	IMPLICIT NONE
60	INTEGER :: ind
61
62	ALLOCATE(vertex_glo(1-halo:iim_glo+halo,1-halo:jjm_glo+halo,nb_face))
63	ncell_glo=10(iim_glo-2)(jjm_glo-2)+(iim_glo-2)*20+12
64	ALLOCATE(cell_glo(ncell_glo))
65	ALLOCATE(tab_index(nb_face,nb_face,0:5))
66	ALLOCATE(edge_glo(ncell_glo*3))
67
68	DO ind=1,ncell_glo
69	cell_glo(ind)%neighbour(:)=0
70	ENDDO
71
72	END SUBROUTINE allocate_metric
73
74	SUBROUTINE Compute_face
75	USE spherical_geom_mod
76	IMPLICIT NONE
77
78	REAL(rstd) :: len_edge
79	REAL(rstd) :: rot=0.
80	INTEGER :: nf,i,j
81	REAL(rstd),DIMENSION(3) :: p1,p2,p3
82	REAL(rstd) :: d1,d2,d3
83
84
85	len_edge=acos(cos(Pi/5)cos(2Pi/5)/(sin(Pi/5)sin(2Pi/5)))
86
87	DO nf=1,nb_face/2
88	CALL lonlat2xyz(0.,Pi/2,vertex_glo(1,1,nf)%xyz)
89	CALL lonlat2xyz(rot+(nf-1)2Pi/5,Pi/2-len_edge,vertex_glo(iim_glo,1,nf)%xyz)
90	CALL lonlat2xyz(rot+2Pi/5+(nf-1)2*Pi/5,Pi/2-len_edge,vertex_glo(1,jjm_glo,nf)%xyz)
91	CALL lonlat2xyz(rot+Pi/5+(nf-1)2Pi/5,-Pi/2+len_edge,vertex_glo(iim_glo,jjm_glo,nf)%xyz)
92
93	CALL lonlat2xyz(0.,-Pi/2,vertex_glo(1,1,nf+nb_face/2)%xyz)
94	CALL lonlat2xyz(rot+Pi/5+(nf-1)2Pi/5,-(Pi/2-len_edge),vertex_glo(1,jjm_glo,nf+nb_face/2)%xyz)
95	CALL lonlat2xyz(rot+Pi/5+2Pi/5+(nf-1)2*Pi/5,-(Pi/2-len_edge),vertex_glo(iim_glo,1,nf+nb_face/2)%xyz)
96	CALL lonlat2xyz(rot+Pi/5+Pi/5+(nf-1)2Pi/5,-(-Pi/2+len_edge),vertex_glo(iim_glo,jjm_glo,nf+nb_face/2)%xyz)
97	ENDDO
98
99	DO nf=1,nb_face
100	DO i=2,iim_glo-1
101	CALL div_arc(vertex_glo(1,1,nf)%xyz,vertex_glo(iim_glo,1,nf)%xyz,(i-1)*1./(iim_glo-1),vertex_glo(i,1,nf)%xyz)
102	CALL div_arc(vertex_glo(1,jjm_glo,nf)%xyz,vertex_glo(iim_glo,jjm_glo,nf)%xyz, &
103	(i-1)*1./(iim_glo-1),vertex_glo(i,jjm_glo,nf)%xyz)
104	ENDDO
105	DO j=2,jjm_glo-1
106	CALL div_arc(vertex_glo(1,1,nf)%xyz,vertex_glo(1,jjm_glo,nf)%xyz,(j-1)*1./(jjm_glo-1),vertex_glo(1,j,nf)%xyz)
107	CALL div_arc(vertex_glo(iim_glo,1,nf)%xyz,vertex_glo(iim_glo,jjm_glo,nf)%xyz,(j-1)*1./(jjm_glo-1), &
108	vertex_glo(iim_glo,j,nf)%xyz)
109	ENDDO
110
111	DO j=2,jjm_glo-1
112	DO i=2,iim_glo-1
113	IF (i+j-1 > iim_glo) THEN
114	CALL div_arc(vertex_glo(iim_glo,i+j-iim_glo,nf)%xyz,vertex_glo(i+j-jjm_glo,jjm_glo,nf)%xyz, &
115	(iim_glo-i)1./(2iim_glo-(i+j)),vertex_glo(i,j,nf)%xyz)
116	ELSE
117	CALL div_arc(vertex_glo(i+j-1,1,nf)%xyz,vertex_glo(1,i+j-1,nf)%xyz,(j-1)*1./(i+j-2),vertex_glo(i,j,nf)%xyz)
118	ENDIF
119	ENDDO
120	ENDDO
121	ENDDO
122
123	CALL dist_cart(vertex_glo(1,1,1)%xyz,vertex_glo(2,1,1)%xyz,d1)
124	CALL dist_cart(vertex_glo(2,1,1)%xyz,vertex_glo(3,1,1)%xyz,d2)
125	CALL dist_cart(vertex_glo(3,1,1)%xyz,vertex_glo(4,1,1)%xyz,d3)
126	CALL div_arc(vertex_glo(1,1,1)%xyz,vertex_glo(3,1,1)%xyz,0.5,p1)
127
128	CALL circumcenter(vertex_glo(1,1,1)%xyz,vertex_glo(2,1,1)%xyz,vertex_glo(1,2,1)%xyz,p1)
129	! CALL Centroide(vertex_glo(1,2,1)%xyz,vertex_glo(2,1,1)%xyz,vertex_glo(1,1,1)%xyz,p1)
130	! CALL Centroide(vertex_glo(1,2,1)%xyz,vertex_glo(1,1,1)%xyz,vertex_glo(2,1,1)%xyz,p1)
131	CALL dist_cart(vertex_glo(1,1,1)%xyz,p1,d1)
132	CALL dist_cart(vertex_glo(2,1,1)%xyz,p1,d2)
133	CALL dist_cart(vertex_glo(1,2,1)%xyz,p1,d3)
134
135	END SUBROUTINE compute_face
136
137	SUBROUTINE Compute_face_projection
138	USE spherical_geom_mod
139	IMPLICIT NONE
140
141	REAL(rstd) :: len_edge
142	REAL(rstd) :: rot=0.
143	INTEGER :: nf,i,j
144	REAL(rstd),DIMENSION(3) :: p1,p2,p3
145	REAL(rstd) :: d1,d2,d3
146
147	len_edge=acos(cos(Pi/5)cos(2Pi/5)/(sin(Pi/5)sin(2Pi/5)))
148
149	DO nf=1,nb_face/2
150	CALL lonlat2xyz(0.,Pi/2,vertex_glo(1,1,nf)%xyz)
151	CALL lonlat2xyz(rot+(nf-1)2Pi/5,Pi/2-len_edge,vertex_glo(iim_glo,1,nf)%xyz)
152	CALL lonlat2xyz(rot+2Pi/5+(nf-1)2*Pi/5,Pi/2-len_edge,vertex_glo(1,jjm_glo,nf)%xyz)
153	CALL lonlat2xyz(rot+Pi/5+(nf-1)2Pi/5,-Pi/2+len_edge,vertex_glo(iim_glo,jjm_glo,nf)%xyz)
154
155	CALL lonlat2xyz(0.,-Pi/2,vertex_glo(1,1,nf+nb_face/2)%xyz)
156	CALL lonlat2xyz(rot+Pi/5+(nf-1)2Pi/5,-(Pi/2-len_edge),vertex_glo(1,jjm_glo,nf+nb_face/2)%xyz)
157	CALL lonlat2xyz(rot+Pi/5+2Pi/5+(nf-1)2*Pi/5,-(Pi/2-len_edge),vertex_glo(iim_glo,1,nf+nb_face/2)%xyz)
158	CALL lonlat2xyz(rot+Pi/5+Pi/5+(nf-1)2Pi/5,-(-Pi/2+len_edge),vertex_glo(iim_glo,jjm_glo,nf+nb_face/2)%xyz)
159	ENDDO
160
161	DO nf=1,nb_face
162	DO i=2,iim_glo-1
163	CALL div_arc_bis(vertex_glo(1,1,nf)%xyz,vertex_glo(iim_glo,1,nf)%xyz,(i-1)*1./(iim_glo-1),vertex_glo(i,1,nf)%xyz)
164	CALL div_arc_bis(vertex_glo(1,jjm_glo,nf)%xyz,vertex_glo(iim_glo,jjm_glo,nf)%xyz, &
165	(i-1)*1./(iim_glo-1),vertex_glo(i,jjm_glo,nf)%xyz)
166	ENDDO
167	DO j=2,jjm_glo-1
168	CALL div_arc_bis(vertex_glo(1,1,nf)%xyz,vertex_glo(1,jjm_glo,nf)%xyz,(j-1)*1./(jjm_glo-1),vertex_glo(1,j,nf)%xyz)
169	CALL div_arc_bis(vertex_glo(iim_glo,1,nf)%xyz,vertex_glo(iim_glo,jjm_glo,nf)%xyz,(j-1)*1./(jjm_glo-1), &
170	vertex_glo(iim_glo,j,nf)%xyz)
171	ENDDO
172
173	DO j=2,jjm_glo-1
174	DO i=2,iim_glo-1
175	IF (i+j-1 > iim_glo) THEN
176	CALL div_arc_bis(vertex_glo(iim_glo,i+j-iim_glo,nf)%xyz,vertex_glo(i+j-jjm_glo,jjm_glo,nf)%xyz, &
177	(iim_glo-i)1./(2iim_glo-(i+j)),vertex_glo(i,j,nf)%xyz)
178	ELSE
179	CALL div_arc_bis(vertex_glo(i+j-1,1,nf)%xyz,vertex_glo(1,i+j-1,nf)%xyz,(j-1)*1./(i+j-2),vertex_glo(i,j,nf)%xyz)
180	ENDIF
181	ENDDO
182	ENDDO
183
184	! DO j=1,jjm_glo
185	! DO i=1,iim_glo
186	! vertex_glo(i,j,nf)%xyz=vertex_glo(i,j,nf)%xyz/sqrt(sum(vertex_glo(i,j,nf)%xyz**2))
187	! ENDDO
188	! ENDDO
189	ENDDO
190
191	! CALL dist_cart(vertex_glo(1,1,1)%xyz,vertex_glo(2,1,1)%xyz,d1)
192	! CALL dist_cart(vertex_glo(2,1,1)%xyz,vertex_glo(3,1,1)%xyz,d2)
193	! CALL dist_cart(vertex_glo(3,1,1)%xyz,vertex_glo(4,1,1)%xyz,d3)
194	! CALL div_arc(vertex_glo(1,1,1)%xyz,vertex_glo(3,1,1)%xyz,0.5,p1)
195	! CALL div_arc(vertex_glo(2,1,1)%xyz,vertex_glo(1,3,1)%xyz,1./3,p1)
196	! CALL div_arc(vertex_glo(1,2,1)%xyz,vertex_glo(3,1,1)%xyz,1./3,p2)
197	! CALL div_arc(vertex_glo(2,2,1)%xyz,vertex_glo(1,1,1)%xyz,1./3,p3)
198	! PRINT *, "dist",d1
199	! PRINT *, "dist",d2
200	! PRINT *, "dist",d3
201	! PRINT *,"dist",vertex_glo(2,1,1)%xyz
202	! PRINT ,"dist",p1/sqrt(sum(p1*2))
203	! CALL Centroide(vertex_glo(1,1,1)%xyz,vertex_glo(2,1,1)%xyz,vertex_glo(1,2,1)%xyz,p1)
204	! CALL Centroide(vertex_glo(1,2,1)%xyz,vertex_glo(2,1,1)%xyz,vertex_glo(1,1,1)%xyz,p1)
205	! CALL Centroide(vertex_glo(1,2,1)%xyz,vertex_glo(1,1,1)%xyz,vertex_glo(2,1,1)%xyz,p1)
206	! CALL dist_cart(vertex_glo(1,1,1)%xyz,p1,d1)
207	! CALL dist_cart(vertex_glo(2,1,1)%xyz,p1,d2)
208	! CALL dist_cart(vertex_glo(1,2,1)%xyz,p1,d3)
209	! PRINT *, "dist",d1
210	! PRINT *, "dist",d2
211	! PRINT *, "dist",d3
212
213	END SUBROUTINE compute_face_projection
214
215
216	SUBROUTINE compute_extended_face_bis
217	IMPLICIT NONE
218
219	INTEGER :: nf
220	INTEGER :: i,j
221	INTEGER :: delta
222	INTEGER :: ind
223	INTEGER :: k,ng
224
225	DO nf=1,nb_face
226	DO i=2,iim_glo-1
227
228	CALL set_neighbour(i,1,nf,ldown-1)
229	CALL set_neighbour(i,1,nf,rdown-1)
230
231	CALL set_neighbour(i,jjm_glo,nf,lup-1)
232	CALL set_neighbour(i,jjm_glo,nf,rup-1)
233
234	ENDDO
235
236	DO j=2,jjm_glo-1
237
238	CALL set_neighbour(1,j,nf,left-1)
239	CALL set_neighbour(1,j,nf,lup-1)
240
241	CALL set_neighbour(iim_glo,j,nf,rdown-1)
242	CALL set_neighbour(iim_glo,j,nf,right-1)
243
244	ENDDO
245
246	CALL set_neighbour(2,0,nf,left-1)
247	CALL set_neighbour(iim_glo,0,nf,right-1)
248	CALL set_neighbour(iim_glo+1,jjm_glo-1,nf,rup-1)
249	CALL set_neighbour(iim_glo-1,jjm_glo+1,nf,right-1)
250	CALL set_neighbour(1,jjm_glo+1,nf,left-1)
251	CALL set_neighbour(0,2,nf,ldown-1)
252
253	CALL set_neighbour(1,0,nf,left-1)
254	CALL set_neighbour(iim_glo,0,nf,right-1)
255	CALL set_neighbour(0,jjm_glo,nf,rup-1)
256	CALL set_neighbour(iim_glo+1,jjm_glo,nf,rup-1)
257	ENDDO
258
259	DO nf=1,nb_face
260	DO j=0,jjm_glo+1
261	DO i=0,iim_glo+1
262	ind=vertex_glo(i,j,nf)%ind
263	delta=vertex_glo(i,j,nf)%delta
264	DO k=0,5
265	ng=MOD(delta+k+6,6)
266	vertex_glo(i,j,nf)%neighbour(k)=cell_glo(ind)%neighbour(ng)
267	ENDDO
268	ENDDO
269	ENDDO
270
271	i=1 ; j=1
272	vertex_glo(i,j,nf)%neighbour(0)=vertex_glo(i+1,j,nf)%ind
273	vertex_glo(i,j,nf)%neighbour(1)=vertex_glo(i,j+1,nf)%ind
274	vertex_glo(i,j,nf)%neighbour(2)=vertex_glo(i-1,j+1,nf)%ind
275	vertex_glo(i,j,nf)%neighbour(3)=vertex_glo(i-1,j,nf)%ind
276	vertex_glo(i,j,nf)%neighbour(4)=vertex_glo(i,j-1,nf)%ind
277	vertex_glo(i,j,nf)%neighbour(5)=vertex_glo(i+1,j-1,nf)%ind
278
279	i=iim_glo ; j=1
280	vertex_glo(i,j,nf)%neighbour(0)=vertex_glo(i+1,j,nf)%ind
281	vertex_glo(i,j,nf)%neighbour(1)=vertex_glo(i,j+1,nf)%ind
282	vertex_glo(i,j,nf)%neighbour(2)=vertex_glo(i-1,j+1,nf)%ind
283	vertex_glo(i,j,nf)%neighbour(3)=vertex_glo(i-1,j,nf)%ind
284	vertex_glo(i,j,nf)%neighbour(4)=vertex_glo(i,j-1,nf)%ind
285	vertex_glo(i,j,nf)%neighbour(5)=vertex_glo(i+1,j-1,nf)%ind
286
287	i=1 ; j=jjm_glo
288	vertex_glo(i,j,nf)%neighbour(0)=vertex_glo(i+1,j,nf)%ind
289	vertex_glo(i,j,nf)%neighbour(1)=vertex_glo(i,j+1,nf)%ind
290	vertex_glo(i,j,nf)%neighbour(2)=vertex_glo(i-1,j+1,nf)%ind
291	vertex_glo(i,j,nf)%neighbour(3)=vertex_glo(i-1,j,nf)%ind
292	vertex_glo(i,j,nf)%neighbour(4)=vertex_glo(i,j-1,nf)%ind
293	vertex_glo(i,j,nf)%neighbour(5)=vertex_glo(i+1,j-1,nf)%ind
294
295	i=iim_glo ; j=jjm_glo
296	vertex_glo(i,j,nf)%neighbour(0)=vertex_glo(i+1,j,nf)%ind
297	vertex_glo(i,j,nf)%neighbour(1)=vertex_glo(i,j+1,nf)%ind
298	vertex_glo(i,j,nf)%neighbour(2)=vertex_glo(i-1,j+1,nf)%ind
299	vertex_glo(i,j,nf)%neighbour(3)=vertex_glo(i-1,j,nf)%ind
300	vertex_glo(i,j,nf)%neighbour(4)=vertex_glo(i,j-1,nf)%ind
301	vertex_glo(i,j,nf)%neighbour(5)=vertex_glo(i+1,j-1,nf)%ind
302	ENDDO
303
304	CONTAINS
305
306	SUBROUTINE set_neighbour(i,j,nf,neighbour)
307	IMPLICIT NONE
308	INTEGER :: i,j,nf,neighbour
309	INTEGER :: in,jn
310	INTEGER :: k,ng
311	INTEGER :: delta
312	INTEGER :: ind,ind2
313
314	SELECT CASE (neighbour)
315	CASE(right-1)
316	in=i+1 ; jn=j
317	CASE(rup-1)
318	in=i ; jn=j+1
319	CASE(lup-1)
320	in=i-1 ; jn=j+1
321	CASE(left-1)
322	in=i-1 ; jn=j
323	CASE(ldown-1)
324	in=i ; jn=j-1
325	CASE(rdown-1)
326	in=i+1; jn=j-1
327	END SELECT
328
329	ind=vertex_glo(i,j,nf)%ind
330	delta=MOD(vertex_glo(i,j,nf)%delta+neighbour+6,6)
331	ind2=cell_glo(ind)%neighbour(delta)
332	DO k=0,5
333	IF (cell_glo(ind2)%neighbour(k)==ind) ng=k
334	ENDDO
335	vertex_glo(in,jn,nf)%ind=ind2
336	vertex_glo(in,jn,nf)%delta=MOD(ng-neighbour+3+6,6)
337	END SUBROUTINE set_neighbour
338
339	END SUBROUTINE compute_extended_face_bis
340
341
342	SUBROUTINE compute_extended_face
343	IMPLICIT NONE
344
345	INTEGER :: nf,nf2
346	INTEGER :: ind,ind2
347	INTEGER :: i,j,h
348
349
350	DO h=1,halo
351	DO nf=1,nb_face
352	DO i=1-h+1,iim_glo+h-1
353	j=1-h+1
354	ind=vertex_glo(i,j,nf)%ind
355	nf2=cell_glo(ind)%assign_face
356	ind2=cell_glo(ind)%neighbour(tab_index(nf,nf2,ldown-1))
357	vertex_glo(i,j-1,nf)%ind=ind2
358	vertex_glo(i,j-1,nf)%xyz=cell_glo(ind2)%xyz
359	ENDDO
360
361	DO j=1-h,jjm_glo+h-1
362	i=iim_glo+h-1
363	ind=vertex_glo(i,j,nf)%ind
364	nf2=cell_glo(ind)%assign_face
365	ind2=cell_glo(ind)%neighbour(tab_index(nf,nf2,right-1))
366	vertex_glo(i+1,j,nf)%ind=ind2
367	vertex_glo(i+1,j,nf)%xyz=cell_glo(ind2)%xyz
368	ENDDO
369
370	DO i=iim_glo+h,1-h+1,-1
371	j=jjm_glo+h-1
372	ind=vertex_glo(i,j,nf)%ind
373	nf2=cell_glo(ind)%assign_face
374	ind2=cell_glo(ind)%neighbour(tab_index(nf,nf2,rup-1))
375	vertex_glo(i,j+1,nf)%ind=ind2
376	vertex_glo(i,j+1,nf)%xyz=cell_glo(ind2)%xyz
377	ENDDO
378
379	DO j=jjm_glo+h,1-h,-1
380	i=1-h+1
381	ind=vertex_glo(i,j,nf)%ind
382	nf2=cell_glo(ind)%assign_face
383	ind2=cell_glo(ind)%neighbour(tab_index(nf,nf2,left-1))
384	vertex_glo(i-1,j,nf)%ind=ind2
385	vertex_glo(i-1,j,nf)%xyz=cell_glo(ind2)%xyz
386	ENDDO
387	ENDDO
388	ENDDO
389
390
391	END SUBROUTINE compute_extended_face
392
393	SUBROUTINE Set_index
394	IMPLICIT NONE
395	INTEGER :: delta
396	INTEGER :: n1,n2,i
397
398	DO n1=1,5
399	DO delta=-2,2
400	DO i=0,5
401	n2=MOD(n1-1-delta+5,5)+1
402	tab_index(n1,n2,i)=MOD(delta+i+6,6)
403	n2=MOD(n1-1+delta+5,5)+1
404	tab_index(n1+5,n2+5,i)=MOD(delta+i+6,6)
405	ENDDO
406	ENDDO
407	ENDDO
408
409	DO n1=1,5
410	DO i=0,5
411	delta=3
412
413	n2=5+n1
414	Tab_index(n1,n2,i)=MOD(delta+i+6,6)
415	Tab_index(n2,n1,i)=MOD(delta+i+6,6)
416
417	n2=5+n1-1
418	IF (n2==5) n2=10
419	Tab_index(n1,n2,i)=MOD(delta+i+6,6)
420	Tab_index(n2,n1,i)=MOD(delta+i+6,6)
421	ENDDO
422	ENDDO
423
424	END SUBROUTINE set_index
425
426
427	SUBROUTINE set_cell
428	IMPLICIT NONE
429	INTEGER :: ind,ind2
430	INTEGER :: nf,nf2,nfm1,nfp1
431	INTEGER :: i,j,k
432	INTEGER :: delta
433
434	ind=0
435
436	!!!!!!!!!!!!!!!!!!!!!!!!!!!!
437	! Association des cellules !
438	!!!!!!!!!!!!!!!!!!!!!!!!!!!!
439
440	! interieur des faces
441	DO nf=1,nb_face
442	DO i=2,iim_glo-1
443	DO j=2,jjm_glo-1
444	ind=ind+1
445	vertex_glo(i,j,nf)%ind=ind
446	cell_glo(ind)%assign_face=nf
447	cell_glo(ind)%assign_i=i
448	cell_glo(ind)%assign_j=j
449	ENDDO
450	ENDDO
451	ENDDO
452
453	! frontieres faces nord
454	DO nf=1,nb_face/2
455
456	nfm1=nf-1
457	IF (nfm1==0) nfm1=nb_face/2
458
459	nf2=nf+nb_face/2
460
461	DO i=2,iim_glo-1
462	ind=ind+1
463	vertex_glo(i,1,nf)%ind=ind
464	vertex_glo(1,i,nfm1)%ind=ind
465	cell_glo(ind)%assign_face=nf
466	cell_glo(ind)%assign_i=i
467	cell_glo(ind)%assign_j=1
468	ENDDO
469
470	DO i=2,iim_glo-1
471	ind=ind+1
472	vertex_glo(i,jjm_glo,nf)%ind=ind
473	vertex_glo(iim_glo-i+1,jjm_glo,nf2)%ind=ind
474
475	cell_glo(ind)%assign_face=nf
476	cell_glo(ind)%assign_i=i
477	cell_glo(ind)%assign_j=jjm_glo
478	ENDDO
479	ENDDO
480
481	! frontieres faces sud
482
483	DO nf=nb_face/2+1,nb_face
484
485	nfp1=nf+1
486	IF (nfp1==nb_face+1) nfp1=nb_face/2+1
487
488	nf2=nf-nb_face/2+1
489	IF (nf2==nb_face/2+1) nf2=1
490
491	DO i=2,iim_glo-1
492	ind=ind+1
493	vertex_glo(i,1,nf)%ind=ind
494	vertex_glo(1,i,nfp1)%ind=ind
495	cell_glo(ind)%assign_face=nf
496	cell_glo(ind)%assign_i=i
497	cell_glo(ind)%assign_j=1
498	ENDDO
499
500	DO j=2,jjm_glo-1
501	ind=ind+1
502	vertex_glo(iim_glo,j,nf)%ind=ind
503	vertex_glo(iim_glo,jjm_glo-j+1,nf2)%ind=ind
504
505	cell_glo(ind)%assign_face=nf
506	cell_glo(ind)%assign_i=iim_glo
507	cell_glo(ind)%assign_j=j
508	ENDDO
509	ENDDO
510
511	! vertex nord (sur 3 faces)
512
513	DO nf=1,nb_face/2
514
515	nfp1=nf+1
516	IF (nfp1==nb_face/2+1) nfp1=1
517
518	nf2=nf+nb_face/2
519
520	ind=ind+1
521	vertex_glo(1,jjm_glo,nf)%ind=ind
522	vertex_glo(iim_glo,1,nfp1)%ind=ind
523	vertex_glo(iim_glo,jjm_glo,nf2)%ind=ind
524	cell_glo(ind)%assign_face=nf
525	cell_glo(ind)%assign_i=1
526	cell_glo(ind)%assign_j=jjm_glo
527	ENDDO
528
529	! vertex sud (sur 3 faces)
530
531	DO nf=nb_face/2+1,nb_face
532
533	nfp1=nf+1
534	IF (nfp1==nb_face+1) nfp1=nb_face/2+1
535
536	nf2=nf-nb_face/2+1
537	IF (nf2==nb_face/2+1) nf2=1
538
539	ind=ind+1
540	vertex_glo(iim_glo,1,nf)%ind=ind
541	vertex_glo(1,jjm_glo,nfp1)%ind=ind
542	vertex_glo(iim_glo,jjm_glo,nf2)%ind=ind
543	cell_glo(ind)%assign_face=nf
544	cell_glo(ind)%assign_i=iim_glo
545	cell_glo(ind)%assign_j=1
546	ENDDO
547
548
549	! vertex nord (sur 5 faces)
550
551	ind=ind+1
552	vertex_glo(1,1,1)%ind=ind
553	vertex_glo(1,1,2)%ind=ind
554	vertex_glo(1,1,3)%ind=ind
555	vertex_glo(1,1,4)%ind=ind
556	vertex_glo(1,1,5)%ind=ind
557	cell_glo(ind)%assign_face=1
558	cell_glo(ind)%assign_i=1
559	cell_glo(ind)%assign_j=1
560
561	! vertex sud (sur 5 faces)
562
563	ind=ind+1
564	vertex_glo(1,1,6)%ind=ind
565	vertex_glo(1,1,7)%ind=ind
566	vertex_glo(1,1,8)%ind=ind
567	vertex_glo(1,1,9)%ind=ind
568	vertex_glo(1,1,10)%ind=ind
569	cell_glo(ind)%assign_face=6
570	cell_glo(ind)%assign_i=1
571	cell_glo(ind)%assign_j=1
572
573	! assignation des coordonnees xyz
574
575	DO ind=1,ncell_glo
576	nf=cell_glo(ind)%assign_face
577	i=cell_glo(ind)%assign_i
578	j=cell_glo(ind)%assign_j
579	cell_glo(ind)%xyz=vertex_glo(i,j,nf)%xyz
580	ENDDO
581
582	!!!!!!!!!!!!!!!!!!!!!!!!!!!
583	! Association des voisins !
584	!!!!!!!!!!!!!!!!!!!!!!!!!!!
585
586	! interieur des faces
587
588	DO nf=1,nb_face
589	DO j=2,jjm_glo-1
590	DO i=2,iim_glo-1
591	ind=vertex_glo(i,j,nf)%ind
592	! right
593	ind2=vertex_glo(i+1,j,nf)%ind
594	nf2=cell_glo(ind2)%assign_face
595	cell_glo(ind)%neighbour(right-1)=ind2
596	cell_glo(ind2)%neighbour(tab_index(nf,nf2,left-1))=ind
597
598	! rup
599	ind2=vertex_glo(i,j+1,nf)%ind
600	nf2=cell_glo(ind2)%assign_face
601	cell_glo(ind)%neighbour(rup-1)=ind2
602	cell_glo(ind2)%neighbour(tab_index(nf,nf2,ldown-1))=ind
603
604	! lup
605	ind2=vertex_glo(i-1,j+1,nf)%ind
606	nf2=cell_glo(ind2)%assign_face
607	cell_glo(ind)%neighbour(lup-1)=ind2
608	cell_glo(ind2)%neighbour(tab_index(nf,nf2,rdown-1))=ind
609
610	! left
611	ind2=vertex_glo(i-1,j,nf)%ind
612	nf2=cell_glo(ind2)%assign_face
613	cell_glo(ind)%neighbour(left-1)=ind2
614	cell_glo(ind2)%neighbour(tab_index(nf,nf2,right-1))=ind
615	! ldown
616	ind2=vertex_glo(i,j-1,nf)%ind
617	nf2=cell_glo(ind2)%assign_face
618	cell_glo(ind)%neighbour(ldown-1)=ind2
619	cell_glo(ind2)%neighbour(tab_index(nf,nf2,rup-1))=ind
620
621	! rdown
622	ind2=vertex_glo(i+1,j-1,nf)%ind
623	nf2=cell_glo(ind2)%assign_face
624	cell_glo(ind)%neighbour(rdown-1)=ind2
625	cell_glo(ind2)%neighbour(tab_index(nf,nf2,lup-1))=ind
626	ENDDO
627	ENDDO
628	ENDDO
629
630	! frontiere face nord
631
632	DO nf=1,nb_face/2
633	DO i=2,iim_glo-1
634	j=1
635	ind=vertex_glo(i,j,nf)%ind
636	! right
637	ind2=vertex_glo(i+1,j,nf)%ind
638	nf2=cell_glo(ind2)%assign_face
639	cell_glo(ind)%neighbour(right-1)=ind2
640	cell_glo(ind2)%neighbour(tab_index(nf,nf2,left-1))=ind
641	! left
642	ind2=vertex_glo(i-1,j,nf)%ind
643	nf2=cell_glo(ind2)%assign_face
644	cell_glo(ind)%neighbour(left-1)=ind2
645	cell_glo(ind2)%neighbour(tab_index(nf,nf2,right-1))=ind
646	! lup
647	ind2=vertex_glo(i-1,j+1,nf)%ind
648	nf2=cell_glo(ind2)%assign_face
649	cell_glo(ind)%neighbour(lup-1)=ind2
650	cell_glo(ind2)%neighbour(tab_index(nf,nf2,rdown-1))=ind
651
652	ENDDO
653
654	DO i=2,iim_glo-1
655	j=jjm_glo
656	ind=vertex_glo(i,j,nf)%ind
657	! right
658	ind2=vertex_glo(i+1,j,nf)%ind
659	nf2=cell_glo(ind2)%assign_face
660	cell_glo(ind)%neighbour(right-1)=ind2
661	cell_glo(ind2)%neighbour(tab_index(nf,nf2,left-1))=ind
662	! left
663	ind2=vertex_glo(i-1,j,nf)%ind
664	nf2=cell_glo(ind2)%assign_face
665	cell_glo(ind)%neighbour(left-1)=ind2
666	cell_glo(ind2)%neighbour(tab_index(nf,nf2,right-1))=ind
667	! rdown
668	ind2=vertex_glo(i+1,j-1,nf)%ind
669	nf2=cell_glo(ind2)%assign_face
670	cell_glo(ind)%neighbour(rdown-1)=ind2
671	cell_glo(ind2)%neighbour(tab_index(nf,nf2,lup-1))=ind
672	ENDDO
673	ENDDO
674
675
676	! frontiere face sud
677
678	DO nf=nb_face/2+1,nb_face
679	DO i=2,iim_glo-1
680	j=1
681	ind=vertex_glo(i,j,nf)%ind
682	! right
683	ind2=vertex_glo(i+1,j,nf)%ind
684	nf2=cell_glo(ind2)%assign_face
685	cell_glo(ind)%neighbour(right-1)=ind2
686	cell_glo(ind2)%neighbour(tab_index(nf,nf2,left-1))=ind
687	! left
688	ind2=vertex_glo(i-1,j,nf)%ind
689	nf2=cell_glo(ind2)%assign_face
690	cell_glo(ind)%neighbour(left-1)=ind2
691	cell_glo(ind2)%neighbour(tab_index(nf,nf2,right-1))=ind
692	! lup
693	ind2=vertex_glo(i-1,j+1,nf)%ind
694	nf2=cell_glo(ind2)%assign_face
695	cell_glo(ind)%neighbour(lup-1)=ind2
696	cell_glo(ind2)%neighbour(tab_index(nf,nf2,rdown-1))=ind
697	ENDDO
698
699	DO j=2,jjm_glo-1
700	i=iim_glo
701	ind=vertex_glo(i,j,nf)%ind
702	! rup
703	ind2=vertex_glo(i,j+1,nf)%ind
704	nf2=cell_glo(ind2)%assign_face
705	cell_glo(ind)%neighbour(rup-1)=ind2
706	cell_glo(ind2)%neighbour(tab_index(nf,nf2,ldown-1))=ind
707	! ldown
708	ind2=vertex_glo(i,j-1,nf)%ind
709	nf2=cell_glo(ind2)%assign_face
710	cell_glo(ind)%neighbour(ldown-1)=ind2
711	cell_glo(ind2)%neighbour(tab_index(nf,nf2,rup-1))=ind
712	! lup
713	ind2=vertex_glo(i-1,j+1,nf)%ind
714	nf2=cell_glo(ind2)%assign_face
715	cell_glo(ind)%neighbour(lup-1)=ind2
716	cell_glo(ind2)%neighbour(tab_index(nf,nf2,rdown-1))=ind
717
718	ENDDO
719	ENDDO
720
721
722	! vertex nord (sur 3 faces)
723
724	DO nf=1,nb_face/2
725	ind=vertex_glo(1,jjm_glo,nf)%ind
726	cell_glo(ind)%neighbour(2)=cell_glo(ind)%neighbour(1)
727	ENDDO
728
729	! vertex sud (sur 3 faces)
730
731	DO nf=nb_face/2+1,nb_face
732
733	ind=vertex_glo(iim_glo,1,nf)%ind
734	cell_glo(ind)%neighbour(5)=cell_glo(ind)%neighbour(0)
735	ENDDO
736
737
738	! vertex nord (sur 5 faces)
739
740	ind=vertex_glo(1,1,1)%ind
741	cell_glo(ind)%neighbour(3)=cell_glo(ind)%neighbour(4)
742
743	! vertex sud (sur 5 faces)
744
745	ind=vertex_glo(1,1,6)%ind
746	cell_glo(ind)%neighbour(3)=cell_glo(ind)%neighbour(2)
747
748
749	!! assignation des delta
750	DO nf=1,nb_face
751	DO j=1,jjm_glo
752	DO i=1,iim_glo
753	ind=vertex_glo(i,j,nf)%ind
754	nf2=cell_glo(ind)%assign_face
755	vertex_glo(i,j,nf)%delta=tab_index(nf,nf2,0)
756	ENDDO
757	ENDDO
758	ENDDO
759
760
761	END SUBROUTINE set_cell
762
763
764	SUBROUTINE set_cell_edge
765	IMPLICIT NONE
766	INTEGER :: ind,k,ind2,ng,ne
767
768	DO ind=1,ncell_glo
769	cell_glo(ind)%edge(:)=0
770	ENDDO
771
772	ne=0
773	DO ind=1,ncell_glo
774	DO k=0,5
775	IF (cell_glo(ind)%edge(k)==0) THEN
776	ind2=cell_glo(ind)%neighbour(k)
777	DO ng=0,5
778	IF (cell_glo(ind2)%neighbour(ng)==ind) THEN
779	IF (cell_glo(ind2)%edge(ng)==0) THEN
780	ne=ne+1
781	cell_glo(ind)%edge(k)=ne
782	edge_glo(ne)%e1=ind
783	edge_glo(ne)%e2=ind2
784	cell_glo(ind2)%edge(ng)=ne
785	ELSE
786	ne=cell_glo(ind2)%edge(ng)
787	cell_glo(ind)%edge(k)=ne
788	ENDIF
789	EXIT
790	ENDIF
791	ENDDO
792	ENDIF
793	ENDDO
794	ENDDO
795
796	END SUBROUTINE set_cell_edge
797
798	SUBROUTINE set_vertex_edge
799	IMPLICIT NONE
800	INTEGER :: nf
801	INTEGER :: i,j
802	INTEGER :: ind,ind2
803	INTEGER :: ne
804	INTEGER :: k,l
805
806	DO nf=1,nb_face
807	DO j=0,jjm_glo+1
808	DO i=0,iim_glo+1
809	ind=vertex_glo(i,j,nf)%ind
810	DO k=0,5
811	ind2=vertex_glo(i,j,nf)%neighbour(k)
812	DO l=0,5
813	ne=cell_glo(ind)%edge(l)
814	IF (edge_glo(ne)%e1==ind .AND. edge_glo(ne)%e2==ind2) THEN
815	vertex_glo(i,j,nf)%ne(k)=1
816	ELSE IF (edge_glo(ne)%e1==ind2 .AND. edge_glo(ne)%e2==ind) THEN
817	vertex_glo(i,j,nf)%ne(k)=-1
818	ENDIF
819	ENDDO
820	ENDDO
821	ENDDO
822	ENDDO
823	ENDDO
824	END SUBROUTINE set_vertex_edge
825
826	SUBROUTINE compute_metric
827	IMPLICIT NONE
828
829	CALL allocate_metric
830	! CALL compute_face
831	CALL compute_face_projection
832	CALL set_index
833	CALL set_cell
834	CALL compute_extended_face_bis
835	CALL set_cell_edge
836	CALL set_vertex_edge
837
838	END SUBROUTINE compute_metric
839
840	END MODULE metric

Note: See TracBrowser for help on using the repository browser.

Download in other formats: