Context Navigation

source: codes/icosagcm/devel/src/diagnostics/wind.F90 @ 612

Last change on this file since 612 was 612, checked in by dubos, 7 years ago
devel : rename directory kernels as kernels_hex
File size: 12.5 KB

Rev	Line
[15]	1	MODULE wind_mod
[584]	2	USE omp_para
	3	USE icosa
	4	IMPLICIT NONE
	5	PRIVATE
[15]	6
[585]	7	PUBLIC :: compute_wind_centered, compute_flux_centered, &
	8	compute_wind_centered_lonlat_compound, compute_wind2d_perp_from_lonlat_compound, &
[584]	9	compute_wind_centered_from_lonlat_compound, compute_wind_perp_from_lonlat_compound, &
[585]	10	un2ulonlat, ulonlat2un
[584]	11
[49]	12	CONTAINS
[15]	13
[151]	14	SUBROUTINE un2ulonlat(f_u, f_ulon, f_ulat)
	15	TYPE(t_field), POINTER :: f_u(:) ! IN : normal velocity components on edges
[584]	16	TYPE(t_field), POINTER :: f_ulon(:), f_ulat(:) ! OUT : velocity reconstructed at hexagons
[151]	17	REAL(rstd),POINTER :: u(:,:), ulon(:,:), ulat(:,:)
	18	INTEGER :: ind
[15]	19
[151]	20	DO ind=1,ndomain
[186]	21	IF (.NOT. assigned_domain(ind)) CYCLE
[151]	22	CALL swap_dimensions(ind)
	23	CALL swap_geometry(ind)
	24	u=f_u(ind)
	25	ulon=f_ulon(ind)
	26	ulat=f_ulat(ind)
	27	CALL compute_un2ulonlat(u,ulon, ulat)
	28	END DO
	29
	30	END SUBROUTINE un2ulonlat
	31
[266]	32	SUBROUTINE ulonlat2un(f_ulon, f_ulat,f_u)
	33	TYPE(t_field), POINTER :: f_ulon(:), f_ulat(:) ! IN : velocity reconstructed at hexagons
	34	TYPE(t_field), POINTER :: f_u(:) ! OUT : normal velocity components on edges
	35
	36	REAL(rstd),POINTER :: u(:,:), ulon(:,:), ulat(:,:)
	37	INTEGER :: ind
	38
	39	DO ind=1,ndomain
	40	IF (.NOT. assigned_domain(ind)) CYCLE
	41	CALL swap_dimensions(ind)
	42	CALL swap_geometry(ind)
	43	u=f_u(ind)
	44	ulon=f_ulon(ind)
	45	ulat=f_ulat(ind)
	46	CALL compute_ulonlat2un(ulon, ulat,u)
	47	END DO
	48
	49	END SUBROUTINE ulonlat2un
[15]	50
	51	SUBROUTINE compute_wind_centered(ue,ucenter)
	52	REAL(rstd) :: ue(3iimjjm,llm)
[584]	53	REAL(rstd) :: ucenter(iim*jjm,llm,3)
[585]	54	INTEGER :: ij,l
[587]	55	REAL(rstd), PARAMETER :: scale=1.
	56	REAL(rstd) :: fac, ue_le, cx,cy,cz, ux,uy,uz
[612]	57	#include "../kernels_hex/wind_centered.k90"
[585]	58	END SUBROUTINE compute_wind_centered
[15]	59
[587]	60	SUBROUTINE compute_flux_centered(scale,ue,ucenter)
	61	REAL(rstd), INTENT(IN) :: scale
[585]	62	REAL(rstd) :: ue(3iimjjm,llm)
	63	REAL(rstd) :: ucenter(iim*jjm,llm,3)
	64	INTEGER :: ij,l
[587]	65	REAL(rstd) :: fac, ue_le, cx,cy,cz, ux,uy,uz
[612]	66	#include "../kernels_hex/flux_centered.k90"
[585]	67	END SUBROUTINE compute_flux_centered
[15]	68
	69
	70	SUBROUTINE compute_wind_on_edge(ue,uedge)
	71	REAL(rstd) :: ue(3iimjjm,llm)
[584]	72	REAL(rstd) :: uedge(3iimjjm,llm,3)
[15]	73
	74	REAL(rstd) :: ut(3iimjjm,llm)
	75	INTEGER :: i,j,ij,l
	76
	77	CALL compute_tangential_compound(ue,ut)
	78
[295]	79	DO l=ll_begin,ll_end
[15]	80	DO j=jj_begin,jj_end
	81	DO i=ii_begin,ii_end
	82	ij=(j-1)*iim+i
[584]	83	uedge(ij+u_right,l,:)=ue(ij+u_right,l)ep_e(ij+u_right,:)ne(ij,right) + ut(ij+u_right,l)et_e(ij+u_right,:)ne(ij,right)
	84	uedge(ij+u_lup,l,:)=ue(ij+u_lup,l)ep_e(ij+u_lup,:)ne(ij,lup) + ut(ij+u_lup,l)et_e(ij+u_lup,:)ne(ij,lup)
	85	uedge(ij+u_ldown,l,:)=ue(ij+u_ldown,l)ep_e(ij+u_ldown,:)ne(ij,ldown) + ut(ij+u_ldown,l)et_e(ij+u_ldown,:)ne(ij,ldown)
[15]	86	ENDDO
	87	ENDDO
	88	ENDDO
	89
	90	END SUBROUTINE compute_wind_on_edge
	91
	92
	93
	94	SUBROUTINE compute_tangential_compound(ue,ut)
	95	REAL(rstd) :: ue(3iimjjm,llm)
	96	REAL(rstd) :: ut(3iimjjm,llm)
	97	INTEGER :: i,j,l,ij
	98
[295]	99	DO l=ll_begin,ll_end
[15]	100	DO j=jj_begin,jj_end
	101	DO i=ii_begin,ii_end
	102	ij=(j-1)*iim+i
	103
	104	ut(ij+u_right,l) = 1/de(ij+u_right) * &
	105	( wee(ij+u_right,1,1)*ue(ij+u_rup,l)+ &
	106	wee(ij+u_right,2,1)*ue(ij+u_lup,l)+ &
	107	wee(ij+u_right,3,1)*ue(ij+u_left,l)+ &
	108	wee(ij+u_right,4,1)*ue(ij+u_ldown,l)+ &
	109	wee(ij+u_right,5,1)*ue(ij+u_rdown,l)+ &
	110	wee(ij+u_right,1,2)*ue(ij+t_right+u_ldown,l)+ &
	111	wee(ij+u_right,2,2)*ue(ij+t_right+u_rdown,l)+ &
	112	wee(ij+u_right,3,2)*ue(ij+t_right+u_right,l)+ &
	113	wee(ij+u_right,4,2)*ue(ij+t_right+u_rup,l)+ &
	114	wee(ij+u_right,5,2)*ue(ij+t_right+u_lup,l) )
	115
	116	ut(ij+u_lup,l) = 1/de(ij+u_lup) * &
	117	( wee(ij+u_lup,1,1)*ue(ij+u_left,l)+ &
	118	wee(ij+u_lup,2,1)*ue(ij+u_ldown,l)+ &
	119	wee(ij+u_lup,3,1)*ue(ij+u_rdown,l)+ &
	120	wee(ij+u_lup,4,1)*ue(ij+u_right,l)+ &
	121	wee(ij+u_lup,5,1)*ue(ij+u_rup,l)+ &
	122	wee(ij+u_lup,1,2)*ue(ij+t_lup+u_right,l)+ &
	123	wee(ij+u_lup,2,2)*ue(ij+t_lup+u_rup,l)+ &
	124	wee(ij+u_lup,3,2)*ue(ij+t_lup+u_lup,l)+ &
	125	wee(ij+u_lup,4,2)*ue(ij+t_lup+u_left,l)+ &
	126	wee(ij+u_lup,5,2)*ue(ij+t_lup+u_ldown,l) )
	127
	128
	129	ut(ij+u_ldown,l) = 1/de(ij+u_ldown) * &
	130	( wee(ij+u_ldown,1,1)*ue(ij+u_rdown,l)+ &
	131	wee(ij+u_ldown,2,1)*ue(ij+u_right,l)+ &
	132	wee(ij+u_ldown,3,1)*ue(ij+u_rup,l)+ &
	133	wee(ij+u_ldown,4,1)*ue(ij+u_lup,l)+ &
	134	wee(ij+u_ldown,5,1)*ue(ij+u_left,l)+ &
	135	wee(ij+u_ldown,1,2)*ue(ij+t_ldown+u_lup,l)+ &
	136	wee(ij+u_ldown,2,2)*ue(ij+t_ldown+u_left,l)+ &
	137	wee(ij+u_ldown,3,2)*ue(ij+t_ldown+u_ldown,l)+ &
	138	wee(ij+u_ldown,4,2)*ue(ij+t_ldown+u_rdown,l)+ &
	139	wee(ij+u_ldown,5,2)*ue(ij+t_ldown+u_right,l) )
	140
	141	ENDDO
	142	ENDDO
	143	ENDDO
	144
	145	END SUBROUTINE compute_tangential_compound
	146
[584]	147	! SUBROUTINE compute_wind_lonlat_compound(u, ulon, ulat)
	148	! REAL(rstd) :: u(3iimjjm,3,llm)
	149	! REAL(rstd) :: ulon(3iimjjm,3,llm)
	150	! REAL(rstd) :: ulat(3iimjjm,3,llm)
	151	!
	152	! INTEGER :: i,j,ij,l
	153	!
	154	!
	155	! DO l=ll_begin,ll_end
	156	! DO j=jj_begin-1,jj_end+1
	157	! DO i=ii_begin-1,ii_end+1
	158	! ij=(j-1)*iim+i
	159	! ulon(ij+u_right,:,l)=sum(u(ij+u_right,:,l)elon_e(ij+u_right,:))elon_e(ij+u_right,:)
	160	! ulon(ij+u_lup,:,l)=sum(u(ij+u_lup,:,l)elon_e(ij+u_lup,:))elon_e(ij+u_lup,:)
	161	! ulon(ij+u_ldown,:,l)=sum(u(ij+u_ldown,:,l)elon_e(ij+u_ldown,:))elon_e(ij+u_ldown,:)
	162	!
	163	! ulat(ij+u_right,:,l)=sum(u(ij+u_right,:,l)elat_e(ij+u_right,:))elat_e(ij+u_right,:)
	164	! ulat(ij+u_lup,:,l)=sum(u(ij+u_lup,:,l)elat_e(ij+u_lup,:))elat_e(ij+u_lup,:)
	165	! ulat(ij+u_ldown,:,l)=sum(u(ij+u_ldown,:,l)elat_e(ij+u_ldown,:))elat_e(ij+u_ldown,:)
	166	!
	167	! ENDDO
	168	! ENDDO
	169	! ENDDO
	170	!
	171	! END SUBROUTINE compute_wind_lonlat_compound
[15]	172
	173	SUBROUTINE compute_wind_from_lonlat_compound(ulon, ulat, u)
[584]	174	REAL(rstd) :: u(3iimjjm,llm,3)
[15]	175	REAL(rstd) :: ulon(3iimjjm,llm)
	176	REAL(rstd) :: ulat(3iimjjm,llm)
	177
	178	INTEGER :: i,j,ij,l
	179
[295]	180	DO l=ll_begin,ll_end
[15]	181	DO j=jj_begin-1,jj_end+1
	182	DO i=ii_begin-1,ii_end+1
	183	ij=(j-1)*iim+i
[584]	184	u(ij+u_right,l,:)=ulon(ij+u_right,l)elon_e(ij+u_right,:)+ ulat(ij+u_right,l)elat_e(ij+u_right,:)
	185	u(ij+u_lup,l,:)=ulon(ij+u_lup,l)elon_e(ij+u_lup,:) + ulat(ij+u_lup,l)elat_e(ij+u_lup,:)
	186	u(ij+u_ldown,l,:)=ulon(ij+u_ldown,l)elon_e(ij+u_ldown,:) + ulat(ij+u_ldown,l)elat_e(ij+u_ldown,:)
[15]	187	ENDDO
	188	ENDDO
	189	ENDDO
	190
	191	END SUBROUTINE compute_wind_from_lonlat_compound
	192
[196]	193	SUBROUTINE compute_wind_centered_from_lonlat_compound(ulon, ulat, u)
[584]	194	REAL(rstd) :: u(iim*jjm,llm,3)
[196]	195	REAL(rstd) :: ulon(iim*jjm,llm)
	196	REAL(rstd) :: ulat(iim*jjm,llm)
	197	INTEGER :: i,j,ij,l
[295]	198	DO l=ll_begin,ll_end
[196]	199	DO j=jj_begin-1,jj_end+1
	200	DO i=ii_begin-1,ii_end+1
	201	ij=(j-1)*iim+i
[584]	202	u(ij,l,:)=ulon(ij,l)elon_i(ij,:)+ ulat(ij,l)elat_i(ij,:)
[196]	203	ENDDO
	204	ENDDO
[584]	205	ENDDO
[196]	206	END SUBROUTINE compute_wind_centered_from_lonlat_compound
	207
[179]	208	SUBROUTINE compute_wind2D_from_lonlat_compound(ulon, ulat, u)
	209	REAL(rstd) :: u(3iimjjm,3)
	210	REAL(rstd) :: ulon(3iimjjm)
	211	REAL(rstd) :: ulat(3iimjjm)
	212
	213	INTEGER :: i,j,ij
	214
	215	DO j=jj_begin-1,jj_end+1
	216	DO i=ii_begin-1,ii_end+1
	217	ij=(j-1)*iim+i
	218	u(ij+u_right,:)=ulon(ij+u_right)elon_e(ij+u_right,:)+ ulat(ij+u_right)elat_e(ij+u_right,:)
	219	u(ij+u_lup,:)=ulon(ij+u_lup)elon_e(ij+u_lup,:) + ulat(ij+u_lup)elat_e(ij+u_lup,:)
	220	u(ij+u_ldown,:)=ulon(ij+u_ldown)elon_e(ij+u_ldown,:) + ulat(ij+u_ldown)elat_e(ij+u_ldown,:)
	221	ENDDO
	222	ENDDO
	223
	224	END SUBROUTINE compute_wind2D_from_lonlat_compound
	225
[36]	226	SUBROUTINE compute_wind_perp_from_lonlat_compound(ulon, ulat, up)
	227	REAL(rstd) :: up(3iimjjm,llm)
	228	REAL(rstd) :: ulon(3iimjjm,llm)
	229	REAL(rstd) :: ulat(3iimjjm,llm)
[584]	230	REAL(rstd) :: u(3iimjjm,llm,3)
[36]	231
	232	INTEGER :: i,j,ij,l
	233
	234	CALL compute_wind_from_lonlat_compound(ulon, ulat, u)
	235
[295]	236	DO l=ll_begin,ll_end
[36]	237	DO j=jj_begin-1,jj_end+1
	238	DO i=ii_begin-1,ii_end+1
	239	ij=(j-1)*iim+i
[584]	240	up(ij+u_right,l)=sum(u(ij+u_right,l,:)*ep_e(ij+u_right,:))
	241	up(ij+u_lup,l)=sum(u(ij+u_lup,l,:)*ep_e(ij+u_lup,:))
	242	up(ij+u_ldown,l)=sum(u(ij+u_ldown,l,:)*ep_e(ij+u_ldown,:))
[36]	243	ENDDO
	244	ENDDO
	245	ENDDO
[15]	246
[36]	247	END SUBROUTINE compute_wind_perp_from_lonlat_compound
	248
[179]	249	SUBROUTINE compute_wind2D_perp_from_lonlat_compound(ulon, ulat, up)
	250	REAL(rstd) :: up(3iimjjm)
	251	REAL(rstd) :: ulon(3iimjjm)
	252	REAL(rstd) :: ulat(3iimjjm)
	253	REAL(rstd) :: u(3iimjjm,3)
	254
	255	INTEGER :: i,j,ij
	256
	257	CALL compute_wind2D_from_lonlat_compound(ulon, ulat, u)
	258	DO j=jj_begin-1,jj_end+1
	259	DO i=ii_begin-1,ii_end+1
	260	ij=(j-1)*iim+i
	261	up(ij+u_right)=sum(u(ij+u_right,:)*ep_e(ij+u_right,:))
	262	up(ij+u_lup)=sum(u(ij+u_lup,:)*ep_e(ij+u_lup,:))
	263	up(ij+u_ldown)=sum(u(ij+u_ldown,:)*ep_e(ij+u_ldown,:))
	264	ENDDO
	265	ENDDO
	266
	267	END SUBROUTINE compute_wind2D_perp_from_lonlat_compound
	268
[15]	269	SUBROUTINE compute_wind_centered_lonlat_compound(uc, ulon, ulat)
[584]	270	REAL(rstd) :: uc(iim*jjm,llm,3)
[38]	271	REAL(rstd) :: ulon(iim*jjm,llm)
	272	REAL(rstd) :: ulat(iim*jjm,llm)
[15]	273
	274	INTEGER :: i,j,ij,l
	275
[295]	276	DO l=ll_begin,ll_end
[15]	277	DO j=jj_begin,jj_end
	278	DO i=ii_begin,ii_end
	279	ij=(j-1)*iim+i
[584]	280	ulon(ij,l)=sum(uc(ij,l,:)*elon_i(ij,:))
	281	ulat(ij,l)=sum(uc(ij,l,:)*elat_i(ij,:))
[15]	282	ENDDO
	283	ENDDO
	284	ENDDO
	285
	286	END SUBROUTINE compute_wind_centered_lonlat_compound
	287
[266]	288	SUBROUTINE compute_wind_centered_from_wind_lonlat_centered(ulon, ulat,uc)
	289	REAL(rstd) :: ulon(iim*jjm,llm)
	290	REAL(rstd) :: ulat(iim*jjm,llm)
[584]	291	REAL(rstd) :: uc(iim*jjm,llm,3)
[266]	292
	293	INTEGER :: i,j,ij,l
	294
	295
[295]	296	DO l=ll_begin,ll_end
[266]	297	DO j=jj_begin,jj_end
	298	DO i=ii_begin,ii_end
	299	ij=(j-1)*iim+i
[584]	300	uc(ij,l,:)=ulon(ij,l)elon_i(ij,:)+ulat(ij,l)elat_i(ij,:)
[266]	301	ENDDO
	302	ENDDO
	303	ENDDO
	304
	305	END SUBROUTINE compute_wind_centered_from_wind_lonlat_centered
	306
[584]	307	SUBROUTINE compute_wind_perp_from_wind_centered(uc,un)
[266]	308
	309	IMPLICIT NONE
[584]	310	REAL(rstd),INTENT(IN) :: uc(iim*jjm,llm,3)
[266]	311	REAL(rstd),INTENT(OUT) :: un(3iimjjm,llm)
	312
	313	INTEGER :: i,j,ij,l
	314
	315
[295]	316	DO l=ll_begin,ll_end
[266]	317	DO j=jj_begin,jj_end
	318	DO i=ii_begin,ii_end
	319	ij=(j-1)*iim+i
[584]	320	un(ij+u_right,l) = sum(0.5(uc(ij,l,:) + uc(ij+t_right,l,:))ep_e(ij+u_right,:))
	321	un(ij+u_lup,l) = sum(0.5(uc(ij,l,:) + uc(ij+t_lup,l,:))ep_e(ij+u_lup,:))
	322	un(ij+u_ldown,l) = sum(0.5(uc(ij,l,:) + uc(ij+t_ldown,l,:))ep_e(ij+u_ldown,:))
[266]	323	ENDDO
	324	ENDDO
	325	ENDDO
	326
	327	END SUBROUTINE compute_wind_perp_from_wind_centered
	328
	329
[151]	330	SUBROUTINE compute_un2ulonlat(un, ulon, ulat)
	331	REAL(rstd),INTENT(IN) :: un(3iimjjm,llm)
	332	REAL(rstd),INTENT(OUT) :: ulon(iim*jjm,llm)
	333	REAL(rstd),INTENT(OUT) :: ulat(iim*jjm,llm)
[15]	334
[584]	335	REAL(rstd) :: uc(iim*jjm,llm,3)
[151]	336
	337	CALL compute_wind_centered(un,uc)
	338	CALL compute_wind_centered_lonlat_compound(uc, ulon, ulat)
	339
	340	END SUBROUTINE compute_un2ulonlat
	341
[266]	342	SUBROUTINE compute_ulonlat2un(ulon, ulat,un)
	343	REAL(rstd),INTENT(IN) :: ulon(iim*jjm,llm)
	344	REAL(rstd),INTENT(IN) :: ulat(iim*jjm,llm)
	345	REAL(rstd),INTENT(OUT) :: un(3iimjjm,llm)
	346
[584]	347	REAL(rstd) :: uc(iim*jjm,llm,3)
[266]	348
	349	CALL compute_wind_centered_from_wind_lonlat_centered(ulon, ulat, uc)
	350	CALL compute_wind_perp_from_wind_centered(uc, un)
	351
	352	END SUBROUTINE compute_ulonlat2un
	353
	354
[15]	355	END MODULE wind_mod

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

Download in other formats: