Context Navigation

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

Last change on this file since 283 was 266, checked in by ymipsl, 10 years ago

Synchronize trunk and Saturn branch.
Merge modification from Saturn branch to trunk

YM

File size: 13.2 KB

Rev	Line
[12]	1	MODULE 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(:,:)
[266]	19	INTEGER,POINTER :: assign_cell_glo(:,:)
[12]	20	INTEGER,POINTER :: assign_i(:,:)
	21	INTEGER,POINTER :: assign_j(:,:)
[21]	22	INTEGER,POINTER :: edge_assign_domain(:,:,:)
	23	INTEGER,POINTER :: edge_assign_i(:,:,:)
	24	INTEGER,POINTER :: edge_assign_j(:,:,:)
	25	INTEGER,POINTER :: edge_assign_pos(:,:,:)
[146]	26	INTEGER,POINTER :: edge_assign_sign(:,:,:)
[12]	27	REAL,POINTER :: xyz(:,:,:)
	28	REAL,POINTER :: neighbour(:,:,:,:)
	29	INTEGER,POINTER :: delta(:,:)
	30	REAL,POINTER :: vertex(:,:,:,:)
	31	LOGICAL,POINTER :: own(:,:)
	32	INTEGER,POINTER :: ne(:,:,:)
	33
	34	INTEGER :: t_right
	35	INTEGER :: t_rup
	36	INTEGER :: t_lup
	37	INTEGER :: t_left
	38	INTEGER :: t_ldown
	39	INTEGER :: t_rdown
	40
	41	INTEGER :: u_right
	42	INTEGER :: u_rup
	43	INTEGER :: u_lup
	44	INTEGER :: u_left
	45	INTEGER :: u_ldown
	46	INTEGER :: u_rdown
	47
	48	INTEGER :: z_rup
	49	INTEGER :: z_up
	50	INTEGER :: z_lup
	51	INTEGER :: z_ldown
	52	INTEGER :: z_down
	53	INTEGER :: z_rdown
	54
	55	INTEGER :: t_pos(6)
	56	INTEGER :: u_pos(6)
	57	INTEGER :: z_pos(6)
	58
	59	END TYPE t_domain
	60
[186]	61	INTEGER,SAVE :: ndomain
	62	INTEGER,SAVE :: ndomain_glo
	63
[12]	64	TYPE(t_domain),SAVE,ALLOCATABLE,TARGET :: domain(:)
[26]	65	TYPE(t_domain),SAVE,ALLOCATABLE,TARGET :: domain_glo(:)
[12]	66
[186]	67	INTEGER,ALLOCATABLE,SAVE :: domglo_rank(:) ! size : ndomain_glo : mpi rank assigned to a domain
	68	INTEGER,ALLOCATABLE,SAVE :: domglo_loc_ind(:) ! size : ndomain_glo : corresponding local indice for a global domain indice
	69	INTEGER,ALLOCATABLE,SAVE :: domloc_glo_ind(:) ! size : ndomain : corresponding global indice for a local domain indice
	70
	71	LOGICAL, ALLOCATABLE, SAVE :: assigned_domain(:) ! size : ndomain : is an omp task is assigned to this domain
	72	!$OMP THREADPRIVATE(assigned_domain)
	73
[12]	74	CONTAINS
	75
	76	SUBROUTINE create_domain
	77	USE grid_param
[151]	78	USE mpipara
[174]	79	USE IOIPSL
[12]	80	IMPLICIT NONE
	81	INTEGER :: ind,nf,ni,nj,i,j
	82	INTEGER :: quotient, rest
[174]	83	INTEGER :: halo_i,halo_j
	84
[12]	85	TYPE(t_domain),POINTER :: d
	86
[26]	87	ndomain_glo=nsplit_insplit_jnb_face
	88	ALLOCATE(domain_glo(ndomain_glo))
	89	ALLOCATE(domglo_rank(ndomain_glo))
	90	ALLOCATE(domglo_loc_ind(ndomain_glo))
[174]	91
	92	halo_i=0
	93	CALL getin("halo_i",halo_i)
	94	halo_j=1
	95	CALL getin("halo_j",halo_j)
[26]	96
[12]	97	ind=0
	98	DO nf=1,nb_face
[21]	99	DO nj=1,nsplit_j
	100	DO ni=1,nsplit_i
[12]	101	ind=ind+1
[26]	102	d=>domain_glo(ind)
[12]	103	quotient=(iim_glo/nsplit_i)
	104	rest=MOD(iim_glo,nsplit_i)
	105	IF (ni-1 < rest) THEN
	106	d%ii_nb=quotient+1
	107	d%ii_begin_glo=(quotient+1)*(ni-1)+1
	108	ELSE
	109	d%ii_nb=quotient
	110	d%ii_begin_glo=(quotient+1)rest+(ni-1-rest) quotient+1
	111	ENDIF
	112	d%ii_end_glo=d%ii_begin_glo+d%ii_nb-1
[21]	113
[151]	114	IF (ni/=1) THEN
[21]	115	d%ii_nb=d%ii_nb+1
[151]	116	d%ii_begin_glo=d%ii_begin_glo-1
[21]	117	ENDIF
	118
[12]	119
	120	quotient=(jjm_glo/nsplit_j)
	121	rest=MOD(jjm_glo,nsplit_j)
	122	IF (nj-1 < rest) THEN
	123	d%jj_nb=quotient+1
	124	d%jj_begin_glo=(quotient+1)*(nj-1)+1
	125	ELSE
	126	d%jj_nb=quotient
	127	d%jj_begin_glo=(quotient+1)rest+(nj-1-rest) quotient+1
	128	ENDIF
[21]	129	d%jj_end_glo=d%jj_begin_glo+d%jj_nb-1
[12]	130
[151]	131	IF (nj/=1) THEN
[21]	132	d%jj_nb=d%jj_nb+1
[151]	133	d%jj_begin_glo=d%jj_begin_glo-1
[21]	134	ENDIF
	135
	136
[174]	137	d%iim=d%ii_nb+2halo + halo_i2
	138	d%jjm=d%jj_nb+2halo + halo_j2
	139	d%ii_begin=halo+1 + halo_i
	140	d%jj_begin=halo+1 + halo_j
[12]	141	d%ii_end=d%ii_begin+d%ii_nb-1
	142	d%jj_end=d%jj_begin+d%jj_nb-1
	143	d%face=nf
	144	ALLOCATE(d%assign_domain(d%iim,d%jjm))
[266]	145	ALLOCATE(d%assign_cell_glo(d%iim,d%jjm))
[12]	146	ALLOCATE(d%assign_i(d%iim,d%jjm))
	147	ALLOCATE(d%assign_j(d%iim,d%jjm))
[21]	148	ALLOCATE(d%edge_assign_domain(0:5,d%iim,d%jjm))
	149	ALLOCATE(d%edge_assign_i(0:5,d%iim,d%jjm))
	150	ALLOCATE(d%edge_assign_j(0:5,d%iim,d%jjm))
	151	ALLOCATE(d%edge_assign_pos(0:5,d%iim,d%jjm))
[146]	152	ALLOCATE(d%edge_assign_sign(0:5,d%iim,d%jjm))
[12]	153	ALLOCATE(d%delta(d%iim,d%jjm))
	154	ALLOCATE(d%xyz(3,d%iim,d%jjm))
	155	ALLOCATE(d%vertex(3,0:5,d%iim,d%jjm))
	156	ALLOCATE(d%neighbour(3,0:5,d%iim,d%jjm))
	157	ALLOCATE(d%own(d%iim,d%jjm))
	158	ALLOCATE(d%ne(0:5,d%iim,d%jjm))
[151]	159
	160	IF (is_mpi_root) PRINT *,"Domain ",ind," : size ",d%ii_nb,"x",d%jj_nb
	161
[12]	162	END DO
	163	END DO
	164	END DO
	165
	166	END SUBROUTINE create_domain
	167
[26]	168	SUBROUTINE copy_domain(d1,d2)
	169	IMPLICIT NONE
	170	INTEGER :: face
	171	TYPE(t_domain),TARGET,INTENT(IN) :: d1
	172	TYPE(t_domain), INTENT(OUT) :: d2
[12]	173
[26]	174	d2%iim = d1%iim
	175	d2%jjm = d1%jjm
	176	d2%ii_begin = d1%ii_begin
	177	d2%jj_begin = d1%jj_begin
	178	d2%ii_end = d1%ii_end
	179	d2%jj_end = d1%jj_end
	180	d2%ii_nb = d1%ii_nb
	181	d2%jj_nb = d1%jj_nb
	182	d2%ii_begin_glo = d1%ii_begin_glo
	183	d2%jj_begin_glo = d1%jj_begin_glo
	184	d2%ii_end_glo = d1%ii_end_glo
	185	d2%jj_end_glo = d1%jj_end_glo
	186	d2%assign_domain => d1%assign_domain
[266]	187	d2%assign_cell_glo => d1%assign_cell_glo
[26]	188	d2%assign_i => d1%assign_i
	189	d2%assign_j => d1%assign_j
	190	d2%edge_assign_domain => d1%edge_assign_domain
	191	d2%edge_assign_i => d1%edge_assign_i
	192	d2%edge_assign_j => d1%edge_assign_j
	193	d2%edge_assign_pos => d1%edge_assign_pos
[146]	194	d2%edge_assign_sign => d1%edge_assign_sign
[26]	195	d2%xyz => d1%xyz
	196	d2%neighbour => d1%neighbour
	197	d2%delta => d1%delta
	198	d2%vertex => d1%vertex
	199	d2%own => d1%own
	200	d2%ne => d1%ne
	201
	202	d2%t_right = d1%t_right
	203	d2%t_rup = d1%t_rup
	204	d2%t_lup = d1%t_lup
	205	d2%t_left = d1%t_left
	206	d2%t_ldown = d1%t_ldown
	207	d2%t_rdown = d1%t_rdown
	208
	209	d2%u_right = d1%u_right
	210	d2%u_rup = d1%u_rup
	211	d2%u_lup = d1%u_lup
	212	d2%u_left = d1%u_left
	213	d2%u_ldown = d1%u_ldown
	214	d2%u_rdown = d1%u_rdown
	215
	216	d2%z_rup = d1%z_rup
	217	d2%z_up = d1%z_up
	218	d2%z_lup = d1%z_lup
	219	d2%z_ldown = d1%z_ldown
	220	d2%z_down = d1%z_down
	221	d2%z_rdown = d1%z_rdown
	222
	223	d2%t_pos = d1%t_pos
	224	d2%u_pos = d1%u_pos
	225	d2%z_pos = d1%z_pos
	226
	227	END SUBROUTINE copy_domain
	228
	229
[12]	230	SUBROUTINE assign_cell
	231	USE metric
	232	IMPLICIT NONE
[21]	233	INTEGER :: ind_d,ind,ind2,e
[146]	234	INTEGER :: nf,nf2
[12]	235	INTEGER :: i,j,k,ii,jj
	236	TYPE(t_domain),POINTER :: d
[21]	237	INTEGER :: delta
[12]	238
	239
[26]	240	DO ind_d=1,ndomain_glo
	241	d=>domain_glo(ind_d)
[12]	242	nf=d%face
	243	DO j=d%jj_begin,d%jj_end
	244	DO i=d%ii_begin,d%ii_end
	245	ii=d%ii_begin_glo-d%ii_begin+i
	246	jj=d%jj_begin_glo-d%jj_begin+j
	247	ind=vertex_glo(ii,jj,nf)%ind
[21]	248	delta=vertex_glo(ii,jj,nf)%delta
[12]	249	IF (cell_glo(ind)%assign_face==nf) THEN
	250	cell_glo(ind)%assign_domain=ind_d
	251	cell_glo(ind)%assign_i=i
	252	cell_glo(ind)%assign_j=j
	253	ENDIF
[21]	254	IF ( i+1 <= d%ii_end ) CALL assign_edge(ind_d,ind,i,j,delta,0)
	255	IF ( j+1 <= d%jj_end ) CALL assign_edge(ind_d,ind,i,j,delta,1)
	256	IF ( i-1 >= d%ii_begin .AND. j+1<=d%jj_end ) CALL assign_edge(ind_d,ind,i,j,delta,2)
	257	IF ( i-1 >= d%ii_begin ) CALL assign_edge(ind_d,ind,i,j,delta,3)
	258	IF ( j-1 >= d%jj_begin ) CALL assign_edge(ind_d,ind,i,j,delta,4)
	259	IF ( i+1 <= d%ii_end .AND. j-1 >=d%jj_begin ) CALL assign_edge(ind_d,ind,i,j,delta,5)
[12]	260	ENDDO
	261	ENDDO
	262	ENDDO
	263
[21]	264
[26]	265	DO ind_d=1,ndomain_glo
	266	d=>domain_glo(ind_d)
[12]	267	nf=d%face
	268	DO j=d%jj_begin-1,d%jj_end+1
	269	DO i=d%ii_begin-1,d%ii_end+1
	270	ii=d%ii_begin_glo-d%ii_begin+i
	271	jj=d%jj_begin_glo-d%jj_begin+j
	272	ind=vertex_glo(ii,jj,nf)%ind
[266]	273	d%assign_cell_glo(i,j) = ind
[12]	274	d%assign_domain(i,j)=cell_glo(ind)%assign_domain
	275	d%assign_i(i,j)=cell_glo(ind)%assign_i
	276	d%assign_j(i,j)=cell_glo(ind)%assign_j
[21]	277	delta=vertex_glo(ii,jj,nf)%delta
[12]	278	d%delta(i,j)=vertex_glo(ii,jj,nf)%delta
	279	DO k=0,5
	280	ind2=vertex_glo(ii,jj,nf)%neighbour(k)
	281	d%neighbour(:,k,i,j)=cell_glo(ind2)%xyz(:)
[146]	282
	283	d%ne(k,i,j)=1-2*MOD(k,2)
	284
[21]	285	e=cell_glo(ind)%edge(MOD(k+delta+6,6))
	286	d%edge_assign_domain(k,i,j)=edge_glo(e)%assign_domain
	287	d%edge_assign_i(k,i,j)=edge_glo(e)%assign_i
	288	d%edge_assign_j(k,i,j)=edge_glo(e)%assign_j
	289	d%edge_assign_pos(k,i,j)=edge_glo(e)%assign_pos
[146]	290	nf2=domain_glo(edge_glo(e)%assign_domain)%face
	291	d%edge_assign_sign(k,i,j)=1-2*MOD(12+tab_index(nf,nf2,0),2)
[149]	292	IF (MOD(6+k+tab_index(nf,nf2,0),6)/=edge_glo(e)%assign_pos .AND. MOD(6+k+tab_index(nf,nf2,0),6) &
	293	/= MOD(edge_glo(e)%assign_pos+3,6)) THEN
[146]	294	d%edge_assign_sign(k,i,j)=-d%edge_assign_sign(k,i,j)
	295	ENDIF
	296
[12]	297	ENDDO
	298	d%xyz(:,i,j)=cell_glo(ind)%xyz(:)
	299	IF (d%assign_domain(i,j)==ind_d) THEN
	300	d%own(i,j)=.TRUE.
	301	ELSE
	302	d%own(i,j)=.FALSE.
	303	ENDIF
	304	ENDDO
	305	ENDDO
[21]	306	ENDDO
	307
	308	CONTAINS
	309
	310	SUBROUTINE assign_edge(ind_d,ind,i,j,delta,k)
	311	INTEGER :: ind_d,ind,i,j,delta,k
	312	INTEGER :: e
	313	e=cell_glo(ind)%edge(MOD(k+delta+6,6))
	314	edge_glo(e)%assign_domain=ind_d
	315	edge_glo(e)%assign_i=i
	316	edge_glo(e)%assign_j=j
	317	edge_glo(e)%assign_pos=k
[146]	318	edge_glo(e)%assign_delta=delta
[151]	319
[21]	320	END SUBROUTINE assign_edge
	321
[12]	322	END SUBROUTINE assign_cell
	323
	324	SUBROUTINE compute_boundary
	325	USE spherical_geom_mod
	326	IMPLICIT NONE
	327	INTEGER :: ind_d
	328	INTEGER :: i,j,k
	329	TYPE(t_domain),POINTER :: d
	330	REAL(rstd) :: ng1(3),ng2(3)
	331
[26]	332	DO ind_d=1,ndomain_glo
	333	d=>domain_glo(ind_d)
[12]	334	DO j=d%jj_begin-1,d%jj_end+1
	335	DO i=d%ii_begin-1,d%ii_end+1
	336	DO k=0,5
	337	ng1=d%neighbour(:,MOD(k,6),i,j)
	338	ng2=d%neighbour(:,MOD(k+1,6),i,j)
	339	IF (sqrt(sum((ng1-ng2)**2))<1e-16) ng2=d%neighbour(:,MOD(k+2,6),i,j)
[15]	340	CALL circumcenter(d%xyz(:,i,j),ng1,ng2,d%vertex(:,k,i,j))
[12]	341	ENDDO
	342	ENDDO
	343	ENDDO
	344	ENDDO
	345	END SUBROUTINE compute_boundary
	346
	347	SUBROUTINE set_neighbour_indice
	348	USE metric
	349	IMPLICIT NONE
	350	INTEGER :: ind_d
	351	TYPE(t_domain),POINTER :: d
	352
[26]	353	DO ind_d=1,ndomain_glo
	354	d=>domain_glo(ind_d)
[12]	355	d%t_right=1
	356	d%t_left=-1
	357	d%t_rup=d%iim
	358	d%t_lup=d%iim-1
	359	d%t_ldown=-d%iim
	360	d%t_rdown=-d%iim+1
	361
	362	d%u_right=0
	363	d%u_lup=d%iim*d%jjm
	364	d%u_ldown=2d%iimd%jjm
	365
	366	d%u_rup=d%t_rup+d%u_ldown
	367	d%u_left=d%t_left+d%u_right
	368	d%u_rdown=d%t_rdown+d%u_lup
	369
	370	d%z_up=0
	371	d%z_down=d%iim*d%jjm
	372	d%z_rup=d%t_rup+d%z_down
	373	d%z_lup=d%t_lup+d%z_down
	374	d%z_ldown=d%t_ldown+d%z_up
	375	d%z_rdown=d%t_rdown+d%z_up
	376
	377	d%t_pos(right)=d%t_right
	378	d%t_pos(rup)=d%t_rup
	379	d%t_pos(lup)=d%t_lup
	380	d%t_pos(left)=d%t_left
	381	d%t_pos(ldown)=d%t_ldown
	382	d%t_pos(rdown)=d%t_rdown
	383
	384	d%u_pos(right)=d%u_right
	385	d%u_pos(rup)=d%u_rup
	386	d%u_pos(lup)=d%u_lup
	387	d%u_pos(left)=d%u_left
	388	d%u_pos(ldown)=d%u_ldown
	389	d%u_pos(rdown)=d%u_rdown
	390
	391	d%z_pos(vrup)=d%z_rup
	392	d%z_pos(vup)=d%z_up
	393	d%z_pos(vlup)=d%z_lup
	394	d%z_pos(vldown)=d%z_ldown
	395	d%z_pos(vdown)=d%z_down
	396	d%z_pos(vrdown)=d%z_rdown
	397
	398	ENDDO
	399
	400	END SUBROUTINE set_neighbour_indice
	401
[26]	402	SUBROUTINE assign_domain
	403	USE mpipara
	404	IMPLICIT NONE
	405	INTEGER :: nb_domain(0:mpi_size-1)
	406	INTEGER :: rank, ind,ind_glo
	407
	408	DO rank=0,mpi_size-1
	409	nb_domain(rank)=ndomain_glo/mpi_size
	410	IF ( rank < MOD(ndomain_glo,mpi_size) ) nb_domain(rank)=nb_domain(rank)+1
	411	ENDDO
	412
	413	ndomain=nb_domain(mpi_rank)
	414	ALLOCATE(domain(ndomain))
	415	ALLOCATE(domloc_glo_ind(ndomain))
	416
	417	rank=0
	418	ind=0
	419	DO ind_glo=1,ndomain_glo
	420	ind=ind+1
	421	domglo_rank(ind_glo)=rank
	422	domglo_loc_ind(ind_glo)=ind
	423	IF (rank==mpi_rank) THEN
	424	CALL copy_domain(domain_glo(ind_glo),domain(ind))
	425	domloc_glo_ind(ind)=ind_glo
	426	ENDIF
[12]	427
[26]	428	IF (ind==nb_domain(rank)) THEN
	429	rank=rank+1
	430	ind=0
	431	ENDIF
	432	ENDDO
[186]	433
	434	!$OMP PARALLEL
	435	CALL assign_domain_omp
	436	!$OMP END PARALLEL
	437
[26]	438
	439	END SUBROUTINE assign_domain
[186]	440
	441	SUBROUTINE assign_domain_omp
	442	USE omp_para
	443	USE mpipara
	444	IMPLICIT NONE
	445	INTEGER :: nb_domain
	446	INTEGER :: rank, ind, i
	447
	448	ALLOCATE(assigned_domain(ndomain))
[26]	449
[186]	450	ind=0
	451	DO rank=0,omp_size-1
	452	nb_domain=ndomain/omp_size
	453	IF ( rank < MOD(ndomain,omp_size) ) nb_domain=nb_domain+1
	454
	455	DO i=1,nb_domain
	456	ind=ind+1
	457	IF (rank==omp_rank) THEN
	458	assigned_domain(ind)=.TRUE.
	459	PRINT *,"Rank ",mpi_rank," task ",rank," local domain : ",ind," global domain ",domloc_glo_ind(ind)
	460	ELSE
	461	assigned_domain(ind)=.FALSE.
	462	ENDIF
	463	ENDDO
	464
	465	ENDDO
	466
	467	END SUBROUTINE assign_domain_omp
	468
	469
[12]	470
	471	SUBROUTINE compute_domain
	472	IMPLICIT NONE
[26]	473	CALL init_domain_param
[12]	474	CALL create_domain
	475	CALL assign_cell
	476	CALL compute_boundary
	477	CALL set_neighbour_indice
[26]	478	CALL assign_domain
[12]	479
	480	END SUBROUTINE compute_domain
	481
	482	END MODULE domain_mod
	483

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

Download in other formats: