Home

Context Navigation

surface_coupler.f90 @ 4882

Last change on this file since 4882 was 4828, checked in by Giersch, 4 years ago
Copyright updated to year 2021, interface pmc_sort removed to accelarate the nesting code
Property svn:keywords set to `Id`
File size: 30.5 KB

Rev	Line
[1682]	1	!> @file surface_coupler.f90
[2000]	2	!------------------------------------------------------------------------------!
[2696]	3	! This file is part of the PALM model system.
[1036]	4	!
[2000]	5	! PALM is free software: you can redistribute it and/or modify it under the
	6	! terms of the GNU General Public License as published by the Free Software
	7	! Foundation, either version 3 of the License, or (at your option) any later
	8	! version.
[1036]	9	!
	10	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
	11	! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
	12	! A PARTICULAR PURPOSE. See the GNU General Public License for more details.
	13	!
	14	! You should have received a copy of the GNU General Public License along with
	15	! PALM. If not, see <http://www.gnu.org/licenses/>.
	16	!
[4828]	17	! Copyright 1997-2021 Leibniz Universitaet Hannover
[2000]	18	!------------------------------------------------------------------------------!
[1036]	19	!
[258]	20	! Current revisions:
[1092]	21	! ------------------
[4671]	22	!
	23	!
[1321]	24	! Former revisions:
	25	! -----------------
	26	! $Id: surface_coupler.f90 4828 2021-01-05 11:21:41Z forkel $
[4671]	27	! Implementation of downward facing USM and LSM surfaces
	28	!
	29	! 4429 2020-02-27 15:24:30Z raasch
[4429]	30	! bugfix: preprocessor directives rearranged for serial mode
	31	!
	32	! 4360 2020-01-07 11:25:50Z suehring
[4182]	33	! Corrected "Former revisions" section
	34	!
	35	! 3655 2019-01-07 16:51:22Z knoop
[3274]	36	! Modularization of all bulk cloud physics code components
[1321]	37	!
[4182]	38	! 109 2007-08-28 15:26:47Z letzel
	39	! Initial revision
	40	!
[102]	41	! Description:
	42	! ------------
[1682]	43	!> Data exchange at the interface between coupled models
[102]	44	!------------------------------------------------------------------------------!
[1682]	45	SUBROUTINE surface_coupler
[4429]	46	#if defined( __parallel )
[1682]	47
[102]	48
[1320]	49	USE arrays_3d, &
[2232]	50	ONLY: pt, rho_ocean, sa, total_2d_a, total_2d_o, u, v
[1320]	51
[3274]	52	USE basic_constants_and_equations_mod, &
	53	ONLY: c_p, l_v
[1427]	54
[1320]	55	USE control_parameters, &
	56	ONLY: coupling_mode, coupling_mode_remote, coupling_topology, &
[2232]	57	humidity, humidity_remote, land_surface, message_string, &
	58	terminate_coupled, terminate_coupled_remote, &
	59	time_since_reference_point, urban_surface
[1320]	60
	61	USE cpulog, &
	62	ONLY: cpu_log, log_point
	63
	64	USE indices, &
	65	ONLY: nbgp, nx, nxl, nxlg, nxr, nxrg, nx_a, nx_o, ny, nyn, nyng, nys, &
	66	nysg, ny_a, ny_o, nzt
	67
	68	USE kinds
	69
[102]	70	USE pegrid
	71
[2232]	72	USE surface_mod, &
	73	ONLY : surf_def_h, surf_lsm_h, surf_type, surf_usm_h
	74
[102]	75	IMPLICIT NONE
	76
[2232]	77	INTEGER(iwp) :: i !< index variable x-direction
	78	INTEGER(iwp) :: j !< index variable y-direction
	79	INTEGER(iwp) :: m !< running index for surface elements
	80
	81	REAL(wp) :: cpw = 4218.0_wp !< heat capacity of water at constant pressure
[1682]	82	REAL(wp) :: time_since_reference_point_rem !<
	83	REAL(wp) :: total_2d(-nbgp:ny+nbgp,-nbgp:nx+nbgp) !<
[102]	84
[2232]	85	REAL(wp), DIMENSION(nysg:nyng,nxlg:nxrg) :: surface_flux !< dummy array for surface fluxes on 2D grid
[1427]	86
[2232]	87
[667]	88	CALL cpu_log( log_point(39), 'surface_coupler', 'start' )
[102]	89
[667]	90
	91
[102]	92	!
[108]	93	!-- In case of model termination initiated by the remote model
	94	!-- (terminate_coupled_remote > 0), initiate termination of the local model.
	95	!-- The rest of the coupler must then be skipped because it would cause an MPI
	96	!-- intercomminucation hang.
	97	!-- If necessary, the coupler will be called at the beginning of the next
	98	!-- restart run.
[667]	99
	100	IF ( coupling_topology == 0 ) THEN
[709]	101	CALL MPI_SENDRECV( terminate_coupled, 1, MPI_INTEGER, target_id, &
	102	0, &
	103	terminate_coupled_remote, 1, MPI_INTEGER, target_id, &
[667]	104	0, comm_inter, status, ierr )
	105	ELSE
	106	IF ( myid == 0) THEN
	107	CALL MPI_SENDRECV( terminate_coupled, 1, MPI_INTEGER, &
	108	target_id, 0, &
	109	terminate_coupled_remote, 1, MPI_INTEGER, &
	110	target_id, 0, &
	111	comm_inter, status, ierr )
	112	ENDIF
[709]	113	CALL MPI_BCAST( terminate_coupled_remote, 1, MPI_INTEGER, 0, comm2d, &
	114	ierr )
[667]	115
	116	ALLOCATE( total_2d_a(-nbgp:ny_a+nbgp,-nbgp:nx_a+nbgp), &
	117	total_2d_o(-nbgp:ny_o+nbgp,-nbgp:nx_o+nbgp) )
	118
	119	ENDIF
	120
[108]	121	IF ( terminate_coupled_remote > 0 ) THEN
[3045]	122	WRITE( message_string, * ) 'remote model "', &
	123	TRIM( coupling_mode_remote ), &
	124	'" terminated', &
[3046]	125	'&with terminate_coupled_remote = ', &
[3045]	126	terminate_coupled_remote, &
[3046]	127	'&local model "', TRIM( coupling_mode ), &
[3045]	128	'" has', &
[3046]	129	'&terminate_coupled = ', &
[667]	130	terminate_coupled
[258]	131	CALL message( 'surface_coupler', 'PA0310', 1, 2, 0, 6, 0 )
[108]	132	RETURN
	133	ENDIF
[667]	134
[291]	135
[108]	136	!
	137	!-- Exchange the current simulated time between the models,
[2232]	138	!-- currently just for total_2d
[709]	139	IF ( coupling_topology == 0 ) THEN
	140
	141	CALL MPI_SEND( time_since_reference_point, 1, MPI_REAL, target_id, 11, &
	142	comm_inter, ierr )
	143	CALL MPI_RECV( time_since_reference_point_rem, 1, MPI_REAL, target_id, &
	144	11, comm_inter, status, ierr )
[667]	145	ELSE
[709]	146
[667]	147	IF ( myid == 0 ) THEN
[709]	148
	149	CALL MPI_SEND( time_since_reference_point, 1, MPI_REAL, target_id, &
	150	11, comm_inter, ierr )
	151	CALL MPI_RECV( time_since_reference_point_rem, 1, MPI_REAL, &
[667]	152	target_id, 11, comm_inter, status, ierr )
[709]	153
[667]	154	ENDIF
[709]	155
	156	CALL MPI_BCAST( time_since_reference_point_rem, 1, MPI_REAL, 0, comm2d, &
	157	ierr )
	158
[667]	159	ENDIF
[102]	160
	161	!
	162	!-- Exchange the interface data
	163	IF ( coupling_mode == 'atmosphere_to_ocean' ) THEN
[667]	164
	165	!
[709]	166	!-- Horizontal grid size and number of processors is equal in ocean and
	167	!-- atmosphere
	168	IF ( coupling_topology == 0 ) THEN
[102]	169
	170	!
[2232]	171	!-- Send heat flux at bottom surface to the ocean. First, transfer from
	172	!-- 1D surface type to 2D grid.
[4671]	173	CALL transfer_1D_to_2D_equal( surf_def_h(0)%shf, surf_lsm_h(0)%shf, &
	174	surf_usm_h(0)%shf )
[2232]	175	CALL MPI_SEND( surface_flux(nysg,nxlg), ngp_xy, MPI_REAL, target_id, &
	176	12, comm_inter, ierr )
[102]	177	!
[2232]	178	!-- Send humidity flux at bottom surface to the ocean. First, transfer
	179	!-- from 1D surface type to 2D grid.
[4671]	180	CALL transfer_1D_to_2D_equal( surf_def_h(0)%qsws, surf_lsm_h(0)%qsws, &
	181	surf_usm_h(0)%qsws )
[667]	182	IF ( humidity ) THEN
[2232]	183	CALL MPI_SEND( surface_flux(nysg,nxlg), ngp_xy, MPI_REAL, &
	184	target_id, 13, comm_inter, ierr )
[667]	185	ENDIF
	186	!
[709]	187	!-- Receive temperature at the bottom surface from the ocean
[2232]	188	CALL MPI_RECV( pt(0,nysg,nxlg), 1, type_xy, target_id, 14, &
[709]	189	comm_inter, status, ierr )
[108]	190	!
[2232]	191	!-- Send the momentum flux (u) at bottom surface to the ocean. First,
	192	!-- transfer from 1D surface type to 2D grid.
[4671]	193	CALL transfer_1D_to_2D_equal( surf_def_h(0)%usws, surf_lsm_h(0)%usws, &
	194	surf_usm_h(0)%usws )
[2232]	195	CALL MPI_SEND( surface_flux(nysg,nxlg), ngp_xy, MPI_REAL, target_id, &
	196	15, comm_inter, ierr )
[102]	197	!
[2232]	198	!-- Send the momentum flux (v) at bottom surface to the ocean. First,
	199	!-- transfer from 1D surface type to 2D grid.
[4671]	200	CALL transfer_1D_to_2D_equal( surf_def_h(0)%vsws, surf_lsm_h(0)%vsws, &
	201	surf_usm_h(0)%vsws )
[2232]	202	CALL MPI_SEND( surface_flux(nysg,nxlg), ngp_xy, MPI_REAL, target_id, &
	203	16, comm_inter, ierr )
[102]	204	!
[709]	205	!-- Receive u at the bottom surface from the ocean
[2232]	206	CALL MPI_RECV( u(0,nysg,nxlg), 1, type_xy, target_id, 17, &
[709]	207	comm_inter, status, ierr )
[667]	208	!
[709]	209	!-- Receive v at the bottom surface from the ocean
[2232]	210	CALL MPI_RECV( v(0,nysg,nxlg), 1, type_xy, target_id, 18, &
[709]	211	comm_inter, status, ierr )
[667]	212	!
	213	!-- Horizontal grid size or number of processors differs between
	214	!-- ocean and atmosphere
	215	ELSE
	216
	217	!
[709]	218	!-- Send heat flux at bottom surface to the ocean
[1353]	219	total_2d_a = 0.0_wp
	220	total_2d = 0.0_wp
[2232]	221	!
	222	!-- Transfer from 1D surface type to 2D grid.
[4671]	223	CALL transfer_1D_to_2D_unequal( surf_def_h(0)%shf, surf_lsm_h(0)%shf, &
	224	surf_usm_h(0)%shf )
[709]	225
[2232]	226	CALL MPI_REDUCE( total_2d, total_2d_a, ngp_a, MPI_REAL, MPI_SUM, 0, &
[709]	227	comm2d, ierr )
	228	CALL interpolate_to_ocean( 12 )
[667]	229	!
[709]	230	!-- Send humidity flux at bottom surface to the ocean
	231	IF ( humidity ) THEN
[1353]	232	total_2d_a = 0.0_wp
	233	total_2d = 0.0_wp
[2232]	234	!
	235	!-- Transfer from 1D surface type to 2D grid.
	236	CALL transfer_1D_to_2D_unequal( surf_def_h(0)%qsws, &
[4671]	237	surf_lsm_h(0)%qsws, &
	238	surf_usm_h(0)%qsws )
[709]	239
	240	CALL MPI_REDUCE( total_2d, total_2d_a, ngp_a, MPI_REAL, MPI_SUM, &
	241	0, comm2d, ierr )
	242	CALL interpolate_to_ocean( 13 )
[667]	243	ENDIF
	244	!
[709]	245	!-- Receive temperature at the bottom surface from the ocean
	246	IF ( myid == 0 ) THEN
[2232]	247	CALL MPI_RECV( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, &
[667]	248	target_id, 14, comm_inter, status, ierr )
	249	ENDIF
	250	CALL MPI_BARRIER( comm2d, ierr )
[709]	251	CALL MPI_BCAST( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, 0, comm2d, &
	252	ierr )
[667]	253	pt(0,nysg:nyng,nxlg:nxrg) = total_2d_a(nysg:nyng,nxlg:nxrg)
	254	!
[709]	255	!-- Send momentum flux (u) at bottom surface to the ocean
[1353]	256	total_2d_a = 0.0_wp
	257	total_2d = 0.0_wp
[2232]	258	!
	259	!-- Transfer from 1D surface type to 2D grid.
[4671]	260	CALL transfer_1D_to_2D_unequal( surf_def_h(0)%usws, surf_lsm_h(0)%usws, &
	261	surf_usm_h(0)%usws )
[709]	262	CALL MPI_REDUCE( total_2d, total_2d_a, ngp_a, MPI_REAL, MPI_SUM, 0, &
	263	comm2d, ierr )
	264	CALL interpolate_to_ocean( 15 )
[667]	265	!
[709]	266	!-- Send momentum flux (v) at bottom surface to the ocean
[1353]	267	total_2d_a = 0.0_wp
	268	total_2d = 0.0_wp
[2232]	269	!
	270	!-- Transfer from 1D surface type to 2D grid.
[4671]	271	CALL transfer_1D_to_2D_unequal( surf_def_h(0)%usws, surf_lsm_h(0)%usws, &
	272	surf_usm_h(0)%usws )
[709]	273	CALL MPI_REDUCE( total_2d, total_2d_a, ngp_a, MPI_REAL, MPI_SUM, 0, &
	274	comm2d, ierr )
	275	CALL interpolate_to_ocean( 16 )
[667]	276	!
[709]	277	!-- Receive u at the bottom surface from the ocean
	278	IF ( myid == 0 ) THEN
[667]	279	CALL MPI_RECV( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, &
[709]	280	target_id, 17, comm_inter, status, ierr )
[667]	281	ENDIF
	282	CALL MPI_BARRIER( comm2d, ierr )
[709]	283	CALL MPI_BCAST( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, 0, comm2d, &
	284	ierr )
[667]	285	u(0,nysg:nyng,nxlg:nxrg) = total_2d_a(nysg:nyng,nxlg:nxrg)
	286	!
[709]	287	!-- Receive v at the bottom surface from the ocean
	288	IF ( myid == 0 ) THEN
[667]	289	CALL MPI_RECV( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, &
[709]	290	target_id, 18, comm_inter, status, ierr )
[667]	291	ENDIF
	292	CALL MPI_BARRIER( comm2d, ierr )
[709]	293	CALL MPI_BCAST( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, 0, comm2d, &
	294	ierr )
[667]	295	v(0,nysg:nyng,nxlg:nxrg) = total_2d_a(nysg:nyng,nxlg:nxrg)
	296
	297	ENDIF
	298
[102]	299	ELSEIF ( coupling_mode == 'ocean_to_atmosphere' ) THEN
	300
	301	!
[667]	302	!-- Horizontal grid size and number of processors is equal
	303	!-- in ocean and atmosphere
	304	IF ( coupling_topology == 0 ) THEN
	305	!
[709]	306	!-- Receive heat flux at the sea surface (top) from the atmosphere
[2232]	307	CALL MPI_RECV( surface_flux(nysg,nxlg), ngp_xy, MPI_REAL, target_id, 12, &
[709]	308	comm_inter, status, ierr )
[2232]	309	CALL transfer_2D_to_1D_equal( surf_def_h(2)%shf )
[102]	310	!
[709]	311	!-- Receive humidity flux from the atmosphere (bottom)
[667]	312	!-- and add it to the heat flux at the sea surface (top)...
	313	IF ( humidity_remote ) THEN
[2232]	314	CALL MPI_RECV( surface_flux(nysg,nxlg), ngp_xy, MPI_REAL, &
[667]	315	target_id, 13, comm_inter, status, ierr )
[2232]	316	CALL transfer_2D_to_1D_equal( surf_def_h(2)%qsws )
[667]	317	ENDIF
	318	!
	319	!-- Send sea surface temperature to the atmosphere model
[709]	320	CALL MPI_SEND( pt(nzt,nysg,nxlg), 1, type_xy, target_id, 14, &
	321	comm_inter, ierr )
[667]	322	!
	323	!-- Receive momentum flux (u) at the sea surface (top) from the atmosphere
[2232]	324	CALL MPI_RECV( surface_flux(nysg,nxlg), ngp_xy, MPI_REAL, target_id, 15, &
[709]	325	comm_inter, status, ierr )
[2232]	326	CALL transfer_2D_to_1D_equal( surf_def_h(2)%usws )
[667]	327	!
	328	!-- Receive momentum flux (v) at the sea surface (top) from the atmosphere
[2232]	329	CALL MPI_RECV( surface_flux(nysg,nxlg), ngp_xy, MPI_REAL, target_id, 16, &
[709]	330	comm_inter, status, ierr )
[2232]	331	CALL transfer_2D_to_1D_equal( surf_def_h(2)%vsws )
[667]	332	!
[709]	333	!-- Send u to the atmosphere
	334	CALL MPI_SEND( u(nzt,nysg,nxlg), 1, type_xy, target_id, 17, &
	335	comm_inter, ierr )
[667]	336	!
[709]	337	!-- Send v to the atmosphere
	338	CALL MPI_SEND( v(nzt,nysg,nxlg), 1, type_xy, target_id, 18, &
	339	comm_inter, ierr )
	340	!
[667]	341	!-- Horizontal gridsize or number of processors differs between
	342	!-- ocean and atmosphere
	343	ELSE
	344	!
[709]	345	!-- Receive heat flux at the sea surface (top) from the atmosphere
	346	IF ( myid == 0 ) THEN
[667]	347	CALL MPI_RECV( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
[709]	348	target_id, 12, comm_inter, status, ierr )
[667]	349	ENDIF
	350	CALL MPI_BARRIER( comm2d, ierr )
[709]	351	CALL MPI_BCAST( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, 0, comm2d, &
	352	ierr )
[2232]	353	CALL transfer_2D_to_1D_unequal( surf_def_h(2)%shf )
[667]	354	!
[709]	355	!-- Receive humidity flux at the sea surface (top) from the atmosphere
	356	IF ( humidity_remote ) THEN
	357	IF ( myid == 0 ) THEN
[667]	358	CALL MPI_RECV( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
[709]	359	target_id, 13, comm_inter, status, ierr )
[667]	360	ENDIF
	361	CALL MPI_BARRIER( comm2d, ierr )
[709]	362	CALL MPI_BCAST( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, 0, &
	363	comm2d, ierr)
[2232]	364	CALL transfer_2D_to_1D_unequal( surf_def_h(2)%qsws )
[667]	365	ENDIF
	366	!
	367	!-- Send surface temperature to atmosphere
[1353]	368	total_2d_o = 0.0_wp
	369	total_2d = 0.0_wp
[667]	370	total_2d(nys:nyn,nxl:nxr) = pt(nzt,nys:nyn,nxl:nxr)
	371
[709]	372	CALL MPI_REDUCE( total_2d, total_2d_o, ngp_o, MPI_REAL, MPI_SUM, 0, &
	373	comm2d, ierr)
	374	CALL interpolate_to_atmos( 14 )
[667]	375	!
[709]	376	!-- Receive momentum flux (u) at the sea surface (top) from the atmosphere
	377	IF ( myid == 0 ) THEN
[667]	378	CALL MPI_RECV( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
[709]	379	target_id, 15, comm_inter, status, ierr )
[667]	380	ENDIF
	381	CALL MPI_BARRIER( comm2d, ierr )
	382	CALL MPI_BCAST( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
[709]	383	0, comm2d, ierr )
[2232]	384	CALL transfer_2D_to_1D_unequal( surf_def_h(2)%usws )
[667]	385	!
[709]	386	!-- Receive momentum flux (v) at the sea surface (top) from the atmosphere
	387	IF ( myid == 0 ) THEN
[667]	388	CALL MPI_RECV( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
[709]	389	target_id, 16, comm_inter, status, ierr )
[667]	390	ENDIF
	391	CALL MPI_BARRIER( comm2d, ierr )
[709]	392	CALL MPI_BCAST( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, 0, comm2d, &
	393	ierr )
[2232]	394	CALL transfer_2D_to_1D_unequal( surf_def_h(2)%vsws )
[667]	395	!
	396	!-- Send u to atmosphere
[1353]	397	total_2d_o = 0.0_wp
	398	total_2d = 0.0_wp
[667]	399	total_2d(nys:nyn,nxl:nxr) = u(nzt,nys:nyn,nxl:nxr)
[709]	400	CALL MPI_REDUCE( total_2d, total_2d_o, ngp_o, MPI_REAL, MPI_SUM, 0, &
	401	comm2d, ierr )
	402	CALL interpolate_to_atmos( 17 )
[667]	403	!
	404	!-- Send v to atmosphere
[1353]	405	total_2d_o = 0.0_wp
	406	total_2d = 0.0_wp
[667]	407	total_2d(nys:nyn,nxl:nxr) = v(nzt,nys:nyn,nxl:nxr)
[709]	408	CALL MPI_REDUCE( total_2d, total_2d_o, ngp_o, MPI_REAL, MPI_SUM, 0, &
	409	comm2d, ierr )
	410	CALL interpolate_to_atmos( 18 )
[667]	411
	412	ENDIF
	413
	414	!
	415	!-- Conversions of fluxes received from atmosphere
	416	IF ( humidity_remote ) THEN
[108]	417	!
[2232]	418	!-- Here top heat flux is still the sum of atmospheric bottom heat fluxes,
[709]	419	!-- * latent heat of vaporization in m2/s2, or 540 cal/g, or 40.65 kJ/mol
	420	!-- /(rho_atm(=1.0)*c_p)
[2232]	421	DO m = 1, surf_def_h(2)%ns
	422	i = surf_def_h(2)%i(m)
	423	j = surf_def_h(2)%j(m)
	424
	425	surf_def_h(2)%shf(m) = surf_def_h(2)%shf(m) + &
[3274]	426	surf_def_h(2)%qsws(m) * l_v / c_p
[709]	427	!
[2232]	428	!-- ...and convert it to a salinity flux at the sea surface (top)
	429	!-- following Steinhorn (1991), JPO 21, pp. 1681-1683:
	430	!-- S'w' = -S * evaporation / ( rho_water * ( 1 - S ) )
	431	surf_def_h(2)%sasws(m) = -1.0_wp * sa(nzt,j,i) * 0.001_wp * &
	432	surf_def_h(2)%qsws(m) / &
	433	( rho_ocean(nzt,j,i) * &
	434	( 1.0_wp - sa(nzt,j,i) * 0.001_wp ) &
	435	)
	436	ENDDO
[108]	437	ENDIF
	438
	439	!
[102]	440	!-- Adjust the kinematic heat flux with respect to ocean density
[2232]	441	!-- (constants are the specific heat capacities for air and water), as well
	442	!-- as momentum fluxes
	443	DO m = 1, surf_def_h(2)%ns
	444	i = surf_def_h(2)%i(m)
	445	j = surf_def_h(2)%j(m)
	446	surf_def_h(2)%shf(m) = surf_def_h(2)%shf(m) / rho_ocean(nzt,j,i) * &
[3274]	447	c_p / cpw
[102]	448
[2232]	449	surf_def_h(2)%usws(m) = surf_def_h(2)%usws(m) / rho_ocean(nzt,j,i)
	450	surf_def_h(2)%vsws(m) = surf_def_h(2)%vsws(m) / rho_ocean(nzt,j,i)
	451	ENDDO
[102]	452
[667]	453	ENDIF
	454
[709]	455	IF ( coupling_topology == 1 ) THEN
[667]	456	DEALLOCATE( total_2d_o, total_2d_a )
	457	ENDIF
	458
	459	CALL cpu_log( log_point(39), 'surface_coupler', 'stop' )
	460
	461
[2232]	462	CONTAINS
	463
	464	! Description:
	465	!------------------------------------------------------------------------------!
	466	!> Data transfer from 1D surface-data type to 2D dummy array for equal
	467	!> grids in atmosphere and ocean.
	468	!------------------------------------------------------------------------------!
	469	SUBROUTINE transfer_1D_to_2D_equal( def_1d, lsm_1d, usm_1d )
	470
	471	IMPLICIT NONE
	472
	473	INTEGER(iwp) :: i !< running index x
	474	INTEGER(iwp) :: j !< running index y
	475	INTEGER(iwp) :: m !< running index surface type
	476
	477	REAL(wp), DIMENSION(1:surf_def_h(0)%ns) :: def_1d !< 1D surface flux, default surfaces
[4671]	478	REAL(wp), DIMENSION(1:surf_lsm_h(0)%ns) :: lsm_1d !< 1D surface flux, natural surfaces
	479	REAL(wp), DIMENSION(1:surf_usm_h(0)%ns) :: usm_1d !< 1D surface flux, urban surfaces
[2232]	480	!
	481	!-- Transfer surface flux at default surfaces to 2D grid
	482	DO m = 1, surf_def_h(0)%ns
	483	i = surf_def_h(0)%i(m)
	484	j = surf_def_h(0)%j(m)
	485	surface_flux(j,i) = def_1d(m)
	486	ENDDO
	487	!
	488	!-- Transfer surface flux at natural surfaces to 2D grid
	489	IF ( land_surface ) THEN
	490	DO m = 1, SIZE(lsm_1d)
[4671]	491	i = surf_lsm_h(0)%i(m)
	492	j = surf_lsm_h(0)%j(m)
[2232]	493	surface_flux(j,i) = lsm_1d(m)
	494	ENDDO
	495	ENDIF
	496	!
	497	!-- Transfer surface flux at natural surfaces to 2D grid
	498	IF ( urban_surface ) THEN
	499	DO m = 1, SIZE(usm_1d)
[4671]	500	i = surf_usm_h(0)%i(m)
	501	j = surf_usm_h(0)%j(m)
[2232]	502	surface_flux(j,i) = usm_1d(m)
	503	ENDDO
	504	ENDIF
	505
	506	END SUBROUTINE transfer_1D_to_2D_equal
	507
	508	! Description:
	509	!------------------------------------------------------------------------------!
	510	!> Data transfer from 2D array for equal grids onto 1D surface-data type
	511	!> array.
	512	!------------------------------------------------------------------------------!
	513	SUBROUTINE transfer_2D_to_1D_equal( def_1d )
	514
	515	IMPLICIT NONE
	516
	517	INTEGER(iwp) :: i !< running index x
	518	INTEGER(iwp) :: j !< running index y
	519	INTEGER(iwp) :: m !< running index surface type
	520
	521	REAL(wp), DIMENSION(1:surf_def_h(2)%ns) :: def_1d !< 1D surface flux, default surfaces
	522	!
	523	!-- Transfer surface flux to 1D surface type, only for default surfaces
	524	DO m = 1, surf_def_h(2)%ns
	525	i = surf_def_h(2)%i(m)
	526	j = surf_def_h(2)%j(m)
	527	def_1d(m) = surface_flux(j,i)
	528	ENDDO
	529
	530	END SUBROUTINE transfer_2D_to_1D_equal
	531
	532	! Description:
	533	!------------------------------------------------------------------------------!
	534	!> Data transfer from 1D surface-data type to 2D dummy array from unequal
	535	!> grids in atmosphere and ocean.
	536	!------------------------------------------------------------------------------!
	537	SUBROUTINE transfer_1D_to_2D_unequal( def_1d, lsm_1d, usm_1d )
	538
	539	IMPLICIT NONE
	540
	541	INTEGER(iwp) :: i !< running index x
	542	INTEGER(iwp) :: j !< running index y
	543	INTEGER(iwp) :: m !< running index surface type
	544
	545	REAL(wp), DIMENSION(1:surf_def_h(0)%ns) :: def_1d !< 1D surface flux, default surfaces
[4671]	546	REAL(wp), DIMENSION(1:surf_lsm_h(0)%ns) :: lsm_1d !< 1D surface flux, natural surfaces
	547	REAL(wp), DIMENSION(1:surf_usm_h(0)%ns) :: usm_1d !< 1D surface flux, urban surfaces
[2232]	548	!
	549	!-- Transfer surface flux at default surfaces to 2D grid. Transfer no
	550	!-- ghost-grid points since total_2d is a global array.
	551	DO m = 1, SIZE(def_1d)
	552	i = surf_def_h(0)%i(m)
	553	j = surf_def_h(0)%j(m)
	554
	555	IF ( i >= nxl .AND. i <= nxr .AND. &
	556	j >= nys .AND. j <= nyn ) THEN
	557	total_2d(j,i) = def_1d(m)
	558	ENDIF
	559	ENDDO
	560	!
	561	!-- Transfer surface flux at natural surfaces to 2D grid
	562	IF ( land_surface ) THEN
	563	DO m = 1, SIZE(lsm_1d)
[4671]	564	i = surf_lsm_h(0)%i(m)
	565	j = surf_lsm_h(0)%j(m)
[2232]	566
	567	IF ( i >= nxl .AND. i <= nxr .AND. &
	568	j >= nys .AND. j <= nyn ) THEN
	569	total_2d(j,i) = lsm_1d(m)
	570	ENDIF
	571	ENDDO
	572	ENDIF
	573	!
	574	!-- Transfer surface flux at natural surfaces to 2D grid
	575	IF ( urban_surface ) THEN
	576	DO m = 1, SIZE(usm_1d)
[4671]	577	i = surf_usm_h(0)%i(m)
	578	j = surf_usm_h(0)%j(m)
[2232]	579
	580	IF ( i >= nxl .AND. i <= nxr .AND. &
	581	j >= nys .AND. j <= nyn ) THEN
	582	total_2d(j,i) = usm_1d(m)
	583	ENDIF
	584	ENDDO
	585	ENDIF
	586
	587	END SUBROUTINE transfer_1D_to_2D_unequal
	588
	589	! Description:
	590	!------------------------------------------------------------------------------!
	591	!> Data transfer from 2D dummy array from unequal grids to 1D surface-data
	592	!> type.
	593	!------------------------------------------------------------------------------!
	594	SUBROUTINE transfer_2D_to_1D_unequal( def_1d )
	595
	596	IMPLICIT NONE
	597
	598	INTEGER(iwp) :: i !< running index x
	599	INTEGER(iwp) :: j !< running index y
	600	INTEGER(iwp) :: m !< running index surface type
	601
	602	REAL(wp), DIMENSION(1:surf_def_h(2)%ns) :: def_1d !< 1D surface flux, default surfaces
	603	!
	604	!-- Transfer 2D surface flux to default surfaces data type. Transfer no
	605	!-- ghost-grid points since total_2d is a global array.
	606	DO m = 1, SIZE(def_1d)
	607	i = surf_def_h(2)%i(m)
	608	j = surf_def_h(2)%j(m)
	609
	610	IF ( i >= nxl .AND. i <= nxr .AND. &
	611	j >= nys .AND. j <= nyn ) THEN
	612	def_1d(m) = total_2d_o(j,i)
	613	ENDIF
	614	ENDDO
	615
	616
	617	END SUBROUTINE transfer_2D_to_1D_unequal
	618
[4429]	619	#endif
[667]	620	END SUBROUTINE surface_coupler
	621
	622
	623
[1682]	624	!------------------------------------------------------------------------------!
	625	! Description:
	626	! ------------
	627	!> @todo Missing subroutine description.
	628	!------------------------------------------------------------------------------!
[4429]	629	#if defined( __parallel )
	630
[709]	631	SUBROUTINE interpolate_to_atmos( tag )
[667]	632
[1320]	633	USE arrays_3d, &
	634	ONLY: total_2d_a, total_2d_o
[667]	635
[1320]	636	USE indices, &
	637	ONLY: nbgp, nx, nx_a, nx_o, ny, ny_a, ny_o
	638
	639	USE kinds
	640
[1324]	641	USE pegrid
[1320]	642
[667]	643	IMPLICIT NONE
	644
[1682]	645	INTEGER(iwp) :: dnx !<
	646	INTEGER(iwp) :: dnx2 !<
	647	INTEGER(iwp) :: dny !<
	648	INTEGER(iwp) :: dny2 !<
	649	INTEGER(iwp) :: i !<
	650	INTEGER(iwp) :: ii !<
	651	INTEGER(iwp) :: j !<
	652	INTEGER(iwp) :: jj !<
[667]	653
[1682]	654	INTEGER(iwp), intent(in) :: tag !<
[1320]	655
[667]	656	CALL MPI_BARRIER( comm2d, ierr )
	657
[709]	658	IF ( myid == 0 ) THEN
	659	!
	660	!-- Cyclic boundary conditions for the total 2D-grid
[667]	661	total_2d_o(-nbgp:-1,:) = total_2d_o(ny+1-nbgp:ny,:)
	662	total_2d_o(:,-nbgp:-1) = total_2d_o(:,nx+1-nbgp:nx)
	663
	664	total_2d_o(ny+1:ny+nbgp,:) = total_2d_o(0:nbgp-1,:)
	665	total_2d_o(:,nx+1:nx+nbgp) = total_2d_o(:,0:nbgp-1)
	666
[102]	667	!
[667]	668	!-- Number of gridpoints of the fine grid within one mesh of the coarse grid
	669	dnx = (nx_o+1) / (nx_a+1)
	670	dny = (ny_o+1) / (ny_a+1)
[102]	671
	672	!
[709]	673	!-- Distance for interpolation around coarse grid points within the fine
	674	!-- grid (note: 2*dnx2 must not be equal with dnx)
[667]	675	dnx2 = 2 * ( dnx / 2 )
	676	dny2 = 2 * ( dny / 2 )
[102]	677
[1353]	678	total_2d_a = 0.0_wp
[102]	679	!
[667]	680	!-- Interpolation from ocean-grid-layer to atmosphere-grid-layer
	681	DO j = 0, ny_a
	682	DO i = 0, nx_a
	683	DO jj = 0, dny2
	684	DO ii = 0, dnx2
	685	total_2d_a(j,i) = total_2d_a(j,i) &
	686	+ total_2d_o(jdny+jj,idnx+ii)
	687	ENDDO
	688	ENDDO
	689	total_2d_a(j,i) = total_2d_a(j,i) / ( ( dnx2 + 1 ) * ( dny2 + 1 ) )
	690	ENDDO
	691	ENDDO
	692	!
[709]	693	!-- Cyclic boundary conditions for atmosphere grid
[667]	694	total_2d_a(-nbgp:-1,:) = total_2d_a(ny_a+1-nbgp:ny_a,:)
	695	total_2d_a(:,-nbgp:-1) = total_2d_a(:,nx_a+1-nbgp:nx_a)
	696
	697	total_2d_a(ny_a+1:ny_a+nbgp,:) = total_2d_a(0:nbgp-1,:)
	698	total_2d_a(:,nx_a+1:nx_a+nbgp) = total_2d_a(:,0:nbgp-1)
	699	!
	700	!-- Transfer of the atmosphere-grid-layer to the atmosphere
[709]	701	CALL MPI_SEND( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, target_id, &
	702	tag, comm_inter, ierr )
[102]	703
	704	ENDIF
	705
[667]	706	CALL MPI_BARRIER( comm2d, ierr )
[102]	707
[4429]	708	END SUBROUTINE interpolate_to_atmos
	709
[880]	710	#endif
	711
[102]	712
[1682]	713	!------------------------------------------------------------------------------!
	714	! Description:
	715	! ------------
	716	!> @todo Missing subroutine description.
	717	!------------------------------------------------------------------------------!
[4429]	718	#if defined( __parallel )
	719
[709]	720	SUBROUTINE interpolate_to_ocean( tag )
[667]	721
[1320]	722	USE arrays_3d, &
	723	ONLY: total_2d_a, total_2d_o
[667]	724
[1320]	725	USE indices, &
	726	ONLY: nbgp, nx, nx_a, nx_o, ny, ny_a, ny_o
	727
	728	USE kinds
	729
[1324]	730	USE pegrid
[1320]	731
[667]	732	IMPLICIT NONE
	733
[1682]	734	INTEGER(iwp) :: dnx !<
	735	INTEGER(iwp) :: dny !<
	736	INTEGER(iwp) :: i !<
	737	INTEGER(iwp) :: ii !<
	738	INTEGER(iwp) :: j !<
	739	INTEGER(iwp) :: jj !<
	740	INTEGER(iwp), intent(in) :: tag !<
[667]	741
[1682]	742	REAL(wp) :: fl !<
	743	REAL(wp) :: fr !<
	744	REAL(wp) :: myl !<
	745	REAL(wp) :: myr !<
[709]	746
[667]	747	CALL MPI_BARRIER( comm2d, ierr )
	748
[709]	749	IF ( myid == 0 ) THEN
[667]	750
	751	!
[709]	752	!-- Number of gridpoints of the fine grid within one mesh of the coarse grid
[667]	753	dnx = ( nx_o + 1 ) / ( nx_a + 1 )
	754	dny = ( ny_o + 1 ) / ( ny_a + 1 )
	755
	756	!
[709]	757	!-- Cyclic boundary conditions for atmosphere grid
[667]	758	total_2d_a(-nbgp:-1,:) = total_2d_a(ny+1-nbgp:ny,:)
	759	total_2d_a(:,-nbgp:-1) = total_2d_a(:,nx+1-nbgp:nx)
	760
	761	total_2d_a(ny+1:ny+nbgp,:) = total_2d_a(0:nbgp-1,:)
	762	total_2d_a(:,nx+1:nx+nbgp) = total_2d_a(:,0:nbgp-1)
	763	!
[709]	764	!-- Bilinear Interpolation from atmosphere grid-layer to ocean grid-layer
[667]	765	DO j = 0, ny
	766	DO i = 0, nx
	767	myl = ( total_2d_a(j+1,i) - total_2d_a(j,i) ) / dny
	768	myr = ( total_2d_a(j+1,i+1) - total_2d_a(j,i+1) ) / dny
	769	DO jj = 0, dny-1
[709]	770	fl = myl*jj + total_2d_a(j,i)
	771	fr = myr*jj + total_2d_a(j,i+1)
[667]	772	DO ii = 0, dnx-1
	773	total_2d_o(jdny+jj,idnx+ii) = ( fr - fl ) / dnx * ii + fl
	774	ENDDO
	775	ENDDO
	776	ENDDO
	777	ENDDO
	778	!
[709]	779	!-- Cyclic boundary conditions for ocean grid
[667]	780	total_2d_o(-nbgp:-1,:) = total_2d_o(ny_o+1-nbgp:ny_o,:)
	781	total_2d_o(:,-nbgp:-1) = total_2d_o(:,nx_o+1-nbgp:nx_o)
	782
	783	total_2d_o(ny_o+1:ny_o+nbgp,:) = total_2d_o(0:nbgp-1,:)
	784	total_2d_o(:,nx_o+1:nx_o+nbgp) = total_2d_o(:,0:nbgp-1)
	785
	786	CALL MPI_SEND( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
	787	target_id, tag, comm_inter, ierr )
	788
	789	ENDIF
	790
	791	CALL MPI_BARRIER( comm2d, ierr )
	792
[4429]	793	END SUBROUTINE interpolate_to_ocean
	794
[880]	795	#endif

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

Context Navigation

source: palm/trunk/SOURCE/surface_coupler.f90 @ 4882

Download in other formats: