Changeset 2606 for palm/trunk/SOURCE/lpm_exchange_horiz.f90
- Timestamp:
- Nov 10, 2017 10:36:31 AM (6 years ago)
- File:
-
- 1 edited
-
palm/trunk/SOURCE/lpm_exchange_horiz.f90 (modified) (20 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/lpm_exchange_horiz.f90
r2330 r2606 25 25 ! ----------------- 26 26 ! $Id$ 27 ! Changed particle box locations: center of particle box now coincides 28 ! with scalar grid point of same index. 29 ! lpm_move_particle now moves inner particles that would leave the core 30 ! to the respective boundary gridbox. Afterwards they get transferred by 31 ! lpm_exchange_horiz. This allows particles to move more than one gridbox 32 ! independent of domain composition. 33 ! Renamed module and subroutines: lpm_pack_arrays_mod -> lpm_pack_and_sort_mod 34 ! lpm_pack_all_arrays -> lpm_sort_in_subboxes, lpm_pack_arrays -> lpm_pack 35 ! lpm_sort -> lpm_sort_timeloop_done 36 ! 37 ! 2330 2017-08-03 14:26:02Z schwenkel 27 38 ! Bugfix: Also for gfortran compilable, function c_sizeof is not used. 28 39 ! … … 128 139 USE kinds 129 140 130 USE lpm_pack_a rrays_mod, &131 ONLY: lpm_pack _arrays141 USE lpm_pack_and_sort_mod, & 142 ONLY: lpm_pack 132 143 133 144 USE netcdf_interface, & … … 255 266 DO n = 1, number_of_particles 256 267 IF ( particles(n)%particle_mask ) THEN 257 i = ( particles(n)%x + 0.5_wp * dx )* ddx268 i = particles(n)%x * ddx 258 269 ! 259 270 !-- Above calculation does not work for indices less than zero 260 IF ( particles(n)%x < -0.5_wp * dx) i = -1271 IF ( particles(n)%x < 0.0_wp) i = -1 261 272 262 273 IF ( i < nxl ) THEN … … 298 309 IF ( particles(n)%particle_mask ) THEN 299 310 300 i = ( particles(n)%x + 0.5_wp * dx )* ddx311 i = particles(n)%x * ddx 301 312 ! 302 313 !-- Above calculation does not work for indices less than zero 303 IF ( particles(n)%x < - 0.5_wp * dx) i = -1314 IF ( particles(n)%x < 0.0_wp ) i = -1 304 315 305 316 IF ( i < nxl ) THEN … … 323 334 deleted_particles = deleted_particles + 1 324 335 325 IF ( trlp(trlp_count)%x >= (nx + 0.5_wp)* dx - 1.0E-12_wp ) THEN336 IF ( trlp(trlp_count)%x >= (nx + 1)* dx - 1.0E-12_wp ) THEN 326 337 trlp(trlp_count)%x = trlp(trlp_count)%x - 1.0E-10_wp 327 338 !++ why is 1 subtracted in next statement??? … … 501 512 DO n = 1, number_of_particles 502 513 IF ( particles(n)%particle_mask ) THEN 503 j = ( particles(n)%y + 0.5_wp * dy )* ddy514 j = particles(n)%y * ddy 504 515 ! 505 516 !-- Above calculation does not work for indices less than zero 506 IF ( particles(n)%y < -0.5_wp * dy) j = -1517 IF ( particles(n)%y < 0.0_wp) j = -1 507 518 508 519 IF ( j < nys ) THEN … … 545 556 IF ( particles(n)%particle_mask ) THEN 546 557 547 j = ( particles(n)%y + 0.5_wp * dy )* ddy558 j = particles(n)%y * ddy 548 559 ! 549 560 !-- Above calculation does not work for indices less than zero 550 IF ( particles(n)%y < -0.5_wp * dy) j = -1561 IF ( particles(n)%y < 0.0_wp ) j = -1 551 562 552 563 IF ( j < nys ) THEN … … 571 582 deleted_particles = deleted_particles + 1 572 583 573 IF ( trsp(trsp_count)%y >= (ny+ 0.5_wp)* dy - 1.0E-12_wp ) THEN584 IF ( trsp(trsp_count)%y >= (ny+1)* dy - 1.0E-12_wp ) THEN 574 585 trsp(trsp_count)%y = trsp(trsp_count)%y - 1.0E-10_wp 575 586 !++ why is 1 subtracted in next statement??? … … 719 730 !-- Apply boundary conditions 720 731 721 IF ( particles(n)%x < -0.5_wp * dx) THEN732 IF ( particles(n)%x < 0.0_wp ) THEN 722 733 723 734 IF ( ibc_par_lr == 0 ) THEN … … 725 736 !-- Cyclic boundary. Relevant coordinate has to be changed. 726 737 particles(n)%x = ( nx + 1 ) * dx + particles(n)%x 727 738 particles(n)%origin_x = ( nx + 1 ) * dx + & 739 particles(n)%origin_x 728 740 ELSEIF ( ibc_par_lr == 1 ) THEN 729 741 ! … … 739 751 ENDIF 740 752 741 ELSEIF ( particles(n)%x >= ( nx + 0.5_wp) * dx ) THEN753 ELSEIF ( particles(n)%x >= ( nx + 1) * dx ) THEN 742 754 743 755 IF ( ibc_par_lr == 0 ) THEN … … 745 757 !-- Cyclic boundary. Relevant coordinate has to be changed. 746 758 particles(n)%x = particles(n)%x - ( nx + 1 ) * dx 759 particles(n)%origin_x = particles(n)%origin_x - & 760 ( nx + 1 ) * dx 747 761 748 762 ELSEIF ( ibc_par_lr == 1 ) THEN … … 773 787 DO n = 1, number_of_particles 774 788 775 IF ( particles(n)%y < -0.5_wp * dy) THEN789 IF ( particles(n)%y < 0.0_wp) THEN 776 790 777 791 IF ( ibc_par_ns == 0 ) THEN … … 779 793 !-- Cyclic boundary. Relevant coordinate has to be changed. 780 794 particles(n)%y = ( ny + 1 ) * dy + particles(n)%y 795 particles(n)%origin_y = ( ny + 1 ) * dy + & 796 particles(n)%origin_y 781 797 782 798 ELSEIF ( ibc_par_ns == 1 ) THEN … … 793 809 ENDIF 794 810 795 ELSEIF ( particles(n)%y >= ( ny + 0.5_wp) * dy ) THEN811 ELSEIF ( particles(n)%y >= ( ny + 1) * dy ) THEN 796 812 797 813 IF ( ibc_par_ns == 0 ) THEN … … 799 815 !-- Cyclic boundary. Relevant coordinate has to be changed. 800 816 particles(n)%y = particles(n)%y - ( ny + 1 ) * dy 817 particles(n)%origin_y = particles(n)%origin_y - & 818 ( ny + 1 ) * dy 801 819 802 820 ELSEIF ( ibc_par_ns == 1 ) THEN … … 866 884 IF ( .NOT. particle_array(n)%particle_mask ) CYCLE 867 885 868 ip = ( particle_array(n)%x + 0.5_wp * dx )* ddx869 jp = ( particle_array(n)%y + 0.5_wp * dy )* ddy870 kp = 886 ip = particle_array(n)%x * ddx 887 jp = particle_array(n)%y * ddy 888 kp = particle_array(n)%z / dz + 1 + offset_ocean_nzt 871 889 872 890 IF ( ip >= nxl .AND. ip <= nxr .AND. jp >= nys .AND. jp <= nyn & … … 880 898 CALL realloc_particles_array (ip,jp,kp) 881 899 ELSE 882 CALL lpm_pack _arrays900 CALL lpm_pack 883 901 prt_count(kp,jp,ip) = number_of_particles 884 902 pindex = prt_count(kp,jp,ip)+1 … … 1014 1032 INTEGER(iwp) :: kp !< index variable along z 1015 1033 INTEGER(iwp) :: n !< index variable for particle array 1016 INTEGER(iwp) :: np_old_cell !< number of particles per grid box before moving 1017 INTEGER(iwp) :: n_start !< start index 1034 INTEGER(iwp) :: np_before_move !< number of particles per grid box before moving 1018 1035 INTEGER(iwp) :: pindex !< dummy argument for number of new particle per grid box 1019 1036 1020 LOGICAL :: pack_done !< 1021 1022 TYPE(particle_type), DIMENSION(:), POINTER :: particles_old_cell !< particles before moving 1037 TYPE(particle_type), DIMENSION(:), POINTER :: particles_before_move !< particles before moving 1023 1038 1024 1039 CALL cpu_log( log_point_s(41), 'lpm_move_particle', 'start' ) … … 1028 1043 DO kp = nzb+1, nzt 1029 1044 1030 np_old_cell = prt_count(kp,jp,ip) 1031 IF ( np_old_cell <= 0 ) CYCLE 1032 particles_old_cell => grid_particles(kp,jp,ip)%particles(1:np_old_cell) 1033 n_start = -1 1045 np_before_move = prt_count(kp,jp,ip) 1046 IF ( np_before_move <= 0 ) CYCLE 1047 particles_before_move => grid_particles(kp,jp,ip)%particles(1:np_before_move) 1034 1048 1035 DO n = 1, np_old_cell 1036 i = ( particles_old_cell(n)%x + 0.5_wp * dx ) * ddx 1037 j = ( particles_old_cell(n)%y + 0.5_wp * dy ) * ddy 1038 k = particles_old_cell(n)%z / dz + 1 + offset_ocean_nzt 1049 DO n = 1, np_before_move 1050 i = particles_before_move(n)%x * ddx 1051 j = particles_before_move(n)%y * ddy 1052 k = particles_before_move(n)%z / dz + 1 + offset_ocean_nzt 1053 1054 !-- For lpm_exchange_horiz to work properly particles need to be moved to the outermost gridboxes 1055 !-- of the respective processor. If the particle index is inside the processor the following lines 1056 !-- will not change the index 1057 i = MIN ( i , nxr ) 1058 i = MAX ( i , nxl ) 1059 j = MIN ( j , nyn ) 1060 j = MAX ( j , nys ) 1061 k = MIN ( k , nzt ) 1062 k = MAX ( k , nzb+1 ) 1039 1063 ! 1040 1064 !-- Check, if particle has moved to another grid cell. 1041 1065 IF ( i /= ip .OR. j /= jp .OR. k /= kp ) THEN 1042 ! 1043 !-- The particle has moved to another grid cell. Now check, if 1044 !-- particle stays on the same processor. 1045 IF ( i >= nxl .AND. i <= nxr .AND. j >= nys .AND. & 1046 j <= nyn .AND. k >= nzb+1 .AND. k <= nzt) THEN 1047 ! 1048 !-- If the particle stays on the same processor, the particle 1049 !-- will be added to the particle array of the new processor. 1050 number_of_particles = prt_count(k,j,i) 1051 particles => grid_particles(k,j,i)%particles(1:number_of_particles) 1052 1053 pindex = prt_count(k,j,i)+1 1054 IF ( pindex > SIZE(grid_particles(k,j,i)%particles) ) & 1055 THEN 1056 n_start = n 1057 EXIT 1058 ENDIF 1059 1060 grid_particles(k,j,i)%particles(pindex) = particles_old_cell(n) 1061 prt_count(k,j,i) = pindex 1062 1063 particles_old_cell(n)%particle_mask = .FALSE. 1066 !! 1067 !-- If the particle stays on the same processor, the particle 1068 !-- will be added to the particle array of the new processor. 1069 number_of_particles = prt_count(k,j,i) 1070 particles => grid_particles(k,j,i)%particles(1:number_of_particles) 1071 1072 pindex = prt_count(k,j,i)+1 1073 IF ( pindex > SIZE(grid_particles(k,j,i)%particles) ) & 1074 THEN 1075 CALL realloc_particles_array(i,j,k) 1064 1076 ENDIF 1077 1078 grid_particles(k,j,i)%particles(pindex) = particles_before_move(n) 1079 prt_count(k,j,i) = pindex 1080 1081 particles_before_move(n)%particle_mask = .FALSE. 1065 1082 ENDIF 1066 1083 ENDDO 1067 1084 1068 IF ( n_start >= 0 ) THEN1069 pack_done = .FALSE.1070 DO n = n_start, np_old_cell1071 i = ( particles_old_cell(n)%x + 0.5_wp * dx ) * ddx1072 j = ( particles_old_cell(n)%y + 0.5_wp * dy ) * ddy1073 k = particles_old_cell(n)%z / dz + 1 + offset_ocean_nzt1074 IF ( i /= ip .OR. j /= jp .OR. k /= kp ) THEN1075 !1076 !-- Particle is in different box1077 IF ( i >= nxl .AND. i <= nxr .AND. j >= nys .AND. &1078 j <= nyn .AND. k >= nzb+1 .AND. k <= nzt) THEN1079 !1080 !-- Particle stays on processor1081 number_of_particles = prt_count(k,j,i)1082 particles => grid_particles(k,j,i)%particles(1:number_of_particles)1083 1084 pindex = prt_count(k,j,i)+11085 IF ( pindex > SIZE(grid_particles(k,j,i)%particles) ) &1086 THEN1087 IF ( pack_done ) THEN1088 CALL realloc_particles_array(i,j,k)1089 ELSE1090 CALL lpm_pack_arrays1091 prt_count(k,j,i) = number_of_particles1092 !1093 !-- Set pindex to next index of the modified1094 !-- particle-array (due to lpm_pack_arrays)1095 pindex = prt_count(k,j,i)+11096 !1097 !-- If number of particles in the new grid box1098 !-- exceeds its allocated memory, the particle array1099 !-- will be reallocated1100 IF ( pindex > SIZE(grid_particles(k,j,i)%particles) ) THEN1101 CALL realloc_particles_array(i,j,k)1102 ENDIF1103 1104 pack_done = .TRUE.1105 ENDIF1106 ENDIF1107 1108 grid_particles(k,j,i)%particles(pindex) = particles_old_cell(n)1109 prt_count(k,j,i) = pindex1110 1111 particles_old_cell(n)%particle_mask = .FALSE.1112 ENDIF1113 ENDIF1114 ENDDO1115 ENDIF1116 1085 ENDDO 1117 1086 ENDDO
Note: See TracChangeset
for help on using the changeset viewer.
