[4777] | 1 | ! |
---|
| 2 | ! $Id$ |
---|
| 3 | ! |
---|
| 4 | ! AGRIF (Adaptive Grid Refinement In Fortran) |
---|
| 5 | ! |
---|
| 6 | ! Copyright (C) 2003 Laurent Debreu (Laurent.Debreu@imag.fr) |
---|
| 7 | ! Christophe Vouland (Christophe.Vouland@imag.fr) |
---|
| 8 | ! |
---|
| 9 | ! This program is free software; you can redistribute it and/or modify |
---|
| 10 | ! it under the terms of the GNU General Public License as published by |
---|
| 11 | ! the Free Software Foundation; either version 2 of the License, or |
---|
| 12 | ! (at your option) any later version. |
---|
| 13 | ! |
---|
| 14 | ! This program is distributed in the hope that it will be useful, |
---|
| 15 | ! but WITHOUT ANY WARRANTY; without even the implied warranty of |
---|
| 16 | ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
---|
| 17 | ! GNU General Public License for more details. |
---|
| 18 | ! |
---|
| 19 | ! You should have received a copy of the GNU General Public License |
---|
| 20 | ! along with this program; if not, write to the Free Software |
---|
| 21 | ! Foundation, Inc., 59 Temple Place- Suite 330, Boston, MA 02111-1307, USA. |
---|
| 22 | ! |
---|
| 23 | ! |
---|
| 24 | !> Module Agrif_Save |
---|
| 25 | !! |
---|
| 26 | !! Module for the initialization by copy of the grids created by clustering. |
---|
| 27 | ! |
---|
| 28 | module Agrif_Save |
---|
| 29 | ! |
---|
| 30 | use Agrif_Types |
---|
| 31 | use Agrif_Link |
---|
| 32 | use Agrif_Arrays |
---|
| 33 | use Agrif_Variables |
---|
| 34 | ! |
---|
| 35 | implicit none |
---|
| 36 | ! |
---|
| 37 | contains |
---|
| 38 | ! |
---|
| 39 | !=================================================================================================== |
---|
| 40 | ! subroutine Agrif_deallocate_Arrays |
---|
| 41 | !--------------------------------------------------------------------------------------------------- |
---|
| 42 | subroutine Agrif_deallocate_Arrays ( var ) |
---|
| 43 | !--------------------------------------------------------------------------------------------------- |
---|
| 44 | type(Agrif_Variable), pointer :: var |
---|
| 45 | ! |
---|
| 46 | if (allocated(var%array1)) deallocate(var%array1) |
---|
| 47 | if (allocated(var%array2)) deallocate(var%array2) |
---|
| 48 | if (allocated(var%array3)) deallocate(var%array3) |
---|
| 49 | if (allocated(var%array4)) deallocate(var%array4) |
---|
| 50 | if (allocated(var%array5)) deallocate(var%array5) |
---|
| 51 | if (allocated(var%array6)) deallocate(var%array6) |
---|
| 52 | ! |
---|
| 53 | if (allocated(var%darray1)) deallocate(var%darray1) |
---|
| 54 | if (allocated(var%darray2)) deallocate(var%darray2) |
---|
| 55 | if (allocated(var%darray3)) deallocate(var%darray3) |
---|
| 56 | if (allocated(var%darray4)) deallocate(var%darray4) |
---|
| 57 | if (allocated(var%darray5)) deallocate(var%darray5) |
---|
| 58 | if (allocated(var%darray6)) deallocate(var%darray6) |
---|
| 59 | ! |
---|
| 60 | if (allocated(var%sarray1)) deallocate(var%sarray1) |
---|
| 61 | if (allocated(var%sarray2)) deallocate(var%sarray2) |
---|
| 62 | if (allocated(var%sarray3)) deallocate(var%sarray3) |
---|
| 63 | if (allocated(var%sarray4)) deallocate(var%sarray4) |
---|
| 64 | if (allocated(var%sarray5)) deallocate(var%sarray5) |
---|
| 65 | if (allocated(var%sarray6)) deallocate(var%sarray6) |
---|
| 66 | ! |
---|
| 67 | ! |
---|
| 68 | if (associated(var%oldvalues2D)) deallocate(var%oldvalues2D) |
---|
| 69 | if (allocated(var%posvar)) deallocate(var%posvar) |
---|
| 70 | if (allocated(var%interptab)) deallocate(var%interptab) |
---|
| 71 | if (allocated(var%coords)) deallocate(var%coords) |
---|
| 72 | !--------------------------------------------------------------------------------------------------- |
---|
| 73 | end subroutine Agrif_deallocate_Arrays |
---|
| 74 | !--------------------------------------------------------------------------------------------------- |
---|
| 75 | subroutine Agrif_deallocate_Arrays_c ( var_c ) |
---|
| 76 | !--------------------------------------------------------------------------------------------------- |
---|
| 77 | type(Agrif_Variable_c), pointer :: var_c |
---|
| 78 | ! |
---|
| 79 | if (allocated(var_c%carray1)) deallocate(var_c%carray1) |
---|
| 80 | if (allocated(var_C%carray2)) deallocate(var_c%carray2) |
---|
| 81 | ! |
---|
| 82 | !--------------------------------------------------------------------------------------------------- |
---|
| 83 | end subroutine Agrif_deallocate_Arrays_c |
---|
| 84 | !=================================================================================================== |
---|
| 85 | !--------------------------------------------------------------------------------------------------- |
---|
| 86 | subroutine Agrif_deallocate_Arrays_l ( var_l ) |
---|
| 87 | !--------------------------------------------------------------------------------------------------- |
---|
| 88 | type(Agrif_Variable_l), pointer :: var_l |
---|
| 89 | ! |
---|
| 90 | ! |
---|
| 91 | if (allocated(var_l%larray1)) deallocate(var_l%larray1) |
---|
| 92 | if (allocated(var_l%larray2)) deallocate(var_l%larray2) |
---|
| 93 | if (allocated(var_l%larray3)) deallocate(var_l%larray3) |
---|
| 94 | if (allocated(var_l%larray4)) deallocate(var_l%larray4) |
---|
| 95 | if (allocated(var_l%larray5)) deallocate(var_l%larray5) |
---|
| 96 | if (allocated(var_l%larray6)) deallocate(var_l%larray6) |
---|
| 97 | ! |
---|
| 98 | !--------------------------------------------------------------------------------------------------- |
---|
| 99 | end subroutine Agrif_deallocate_Arrays_l |
---|
| 100 | !--------------------------------------------------------------------------------------------------- |
---|
| 101 | subroutine Agrif_deallocate_Arrays_i ( var_i ) |
---|
| 102 | !--------------------------------------------------------------------------------------------------- |
---|
| 103 | type(Agrif_Variable_i), pointer :: var_i |
---|
| 104 | ! |
---|
| 105 | ! |
---|
| 106 | if (allocated(var_i%iarray1)) deallocate(var_i%iarray1) |
---|
| 107 | if (allocated(var_i%iarray2)) deallocate(var_i%iarray2) |
---|
| 108 | if (allocated(var_i%iarray3)) deallocate(var_i%iarray3) |
---|
| 109 | if (allocated(var_i%iarray4)) deallocate(var_i%iarray4) |
---|
| 110 | if (allocated(var_i%iarray5)) deallocate(var_i%iarray5) |
---|
| 111 | if (allocated(var_i%iarray6)) deallocate(var_i%iarray6) |
---|
| 112 | ! |
---|
| 113 | !--------------------------------------------------------------------------------------------------- |
---|
| 114 | end subroutine Agrif_deallocate_Arrays_i |
---|
| 115 | ! |
---|
| 116 | !=================================================================================================== |
---|
| 117 | ! subroutine Agrif_Free_data_before |
---|
| 118 | !--------------------------------------------------------------------------------------------------- |
---|
| 119 | subroutine Agrif_Free_data_before ( Agrif_Gr ) |
---|
| 120 | !--------------------------------------------------------------------------------------------------- |
---|
| 121 | type(Agrif_Grid), pointer :: Agrif_Gr ! Pointer on the current grid |
---|
| 122 | ! |
---|
| 123 | integer :: i |
---|
| 124 | type(Agrif_Variables_List), pointer :: parcours |
---|
| 125 | type(Agrif_Variable), pointer :: var |
---|
| 126 | type(Agrif_Variable_c), pointer :: var_c |
---|
| 127 | type(Agrif_Variable_r), pointer :: var_r |
---|
| 128 | type(Agrif_Variable_l), pointer :: var_l |
---|
| 129 | type(Agrif_Variable_i), pointer :: var_i |
---|
| 130 | ! |
---|
| 131 | do i = 1,Agrif_NbVariables(0) |
---|
| 132 | ! |
---|
| 133 | var => Agrif_Gr % tabvars(i) |
---|
| 134 | ! |
---|
| 135 | if ( .NOT. Agrif_Mygrid % tabvars(i) % restore ) then |
---|
| 136 | call Agrif_deallocate_Arrays(var) |
---|
| 137 | endif |
---|
| 138 | ! |
---|
| 139 | if (associated(var%list_interp)) then |
---|
| 140 | call Agrif_Free_list_interp(var%list_interp) |
---|
| 141 | endif |
---|
| 142 | ! |
---|
| 143 | enddo |
---|
| 144 | do i=1,Agrif_NbVariables(1) |
---|
| 145 | var_c => Agrif_Gr % tabvars_c(i) |
---|
| 146 | call Agrif_deallocate_Arrays_c(var_c) |
---|
| 147 | enddo |
---|
| 148 | do i=1,Agrif_NbVariables(3) |
---|
| 149 | var_l => Agrif_Gr % tabvars_l(i) |
---|
| 150 | call Agrif_deallocate_Arrays_l(var_l) |
---|
| 151 | enddo |
---|
| 152 | do i=1,Agrif_NbVariables(4) |
---|
| 153 | var_i => Agrif_Gr % tabvars_i(i) |
---|
| 154 | call Agrif_deallocate_Arrays_i(var_i) |
---|
| 155 | enddo |
---|
| 156 | |
---|
| 157 | parcours => Agrif_Gr % variables |
---|
| 158 | |
---|
| 159 | do i = 1,Agrif_Gr%NbVariables |
---|
| 160 | ! |
---|
| 161 | if ( .NOT. parcours%var%root_var % restore ) then |
---|
| 162 | call Agrif_deallocate_Arrays(parcours%var) |
---|
| 163 | endif |
---|
| 164 | ! |
---|
| 165 | if (associated(parcours%var%list_interp)) then |
---|
| 166 | call Agrif_Free_list_interp(parcours%var%list_interp) |
---|
| 167 | endif |
---|
| 168 | ! |
---|
| 169 | if ( .NOT. parcours%var%root_var % restore ) then |
---|
| 170 | deallocate(parcours%var) |
---|
| 171 | endif |
---|
| 172 | ! |
---|
| 173 | parcours => parcours % next |
---|
| 174 | ! |
---|
| 175 | enddo |
---|
| 176 | ! |
---|
| 177 | if ( Agrif_USE_ONLY_FIXED_GRIDS == 0 ) then |
---|
| 178 | if ( Agrif_Probdim == 1 ) deallocate(Agrif_Gr%tabpoint1D) |
---|
| 179 | if ( Agrif_Probdim == 2 ) deallocate(Agrif_Gr%tabpoint2D) |
---|
| 180 | if ( Agrif_Probdim == 3 ) deallocate(Agrif_Gr%tabpoint3D) |
---|
| 181 | endif |
---|
| 182 | !--------------------------------------------------------------------------------------------------- |
---|
| 183 | end subroutine Agrif_Free_data_before |
---|
| 184 | !=================================================================================================== |
---|
| 185 | ! |
---|
| 186 | !=================================================================================================== |
---|
| 187 | ! subroutine Agrif_Free_list_interp |
---|
| 188 | !--------------------------------------------------------------------------------------------------- |
---|
| 189 | recursive subroutine Agrif_Free_list_interp ( list_interp ) |
---|
| 190 | !--------------------------------------------------------------------------------------------------- |
---|
| 191 | type(Agrif_List_Interp_Loc), pointer :: list_interp |
---|
| 192 | ! |
---|
| 193 | if (associated(list_interp%suiv)) call Agrif_Free_list_interp(list_interp%suiv) |
---|
| 194 | |
---|
| 195 | #if defined AGRIF_MPI |
---|
| 196 | deallocate(list_interp%interp_loc%tab4t) |
---|
| 197 | deallocate(list_interp%interp_loc%memberinall) |
---|
| 198 | deallocate(list_interp%interp_loc%sendtoproc1) |
---|
| 199 | deallocate(list_interp%interp_loc%recvfromproc1) |
---|
| 200 | #endif |
---|
| 201 | deallocate(list_interp%interp_loc) |
---|
| 202 | deallocate(list_interp) |
---|
| 203 | nullify(list_interp) |
---|
| 204 | !--------------------------------------------------------------------------------------------------- |
---|
| 205 | end subroutine Agrif_Free_list_interp |
---|
| 206 | !=================================================================================================== |
---|
| 207 | ! |
---|
| 208 | !=================================================================================================== |
---|
| 209 | ! subroutine Agrif_Free_data_after |
---|
| 210 | !--------------------------------------------------------------------------------------------------- |
---|
| 211 | subroutine Agrif_Free_data_after ( Agrif_Gr ) |
---|
| 212 | !--------------------------------------------------------------------------------------------------- |
---|
| 213 | type(Agrif_Grid), pointer :: Agrif_Gr !< Pointer on the current grid |
---|
| 214 | ! |
---|
| 215 | integer :: i |
---|
| 216 | type(Agrif_Variables_List), pointer :: parcours, rootparcours |
---|
| 217 | type(Agrif_Variable), pointer :: var |
---|
| 218 | ! |
---|
| 219 | do i = 1,Agrif_NbVariables(0) |
---|
| 220 | if ( Agrif_Mygrid % tabvars(i) % restore ) then |
---|
| 221 | var => Agrif_Gr%tabvars(i) |
---|
| 222 | call Agrif_deallocate_Arrays(var) |
---|
| 223 | endif |
---|
| 224 | enddo |
---|
| 225 | ! |
---|
| 226 | parcours => Agrif_Gr%variables |
---|
| 227 | rootparcours => Agrif_Mygrid%variables |
---|
| 228 | ! |
---|
| 229 | do i = 1,Agrif_Gr%NbVariables |
---|
| 230 | if (rootparcours%var % restore ) then |
---|
| 231 | call Agrif_deallocate_Arrays(parcours%var) |
---|
| 232 | deallocate(parcours%var) |
---|
| 233 | endif |
---|
| 234 | parcours => parcours % next |
---|
| 235 | rootparcours => rootparcours % next |
---|
| 236 | enddo |
---|
| 237 | !--------------------------------------------------------------------------------------------------- |
---|
| 238 | end subroutine Agrif_Free_data_after |
---|
| 239 | !=================================================================================================== |
---|
| 240 | ! |
---|
| 241 | !=================================================================================================== |
---|
| 242 | ! subroutine Agrif_CopyFromOld_All |
---|
| 243 | ! |
---|
| 244 | !> Called in Agrif_Clustering#Agrif_Init_Hierarchy. |
---|
| 245 | !--------------------------------------------------------------------------------------------------- |
---|
| 246 | recursive subroutine Agrif_CopyFromOld_All ( g, oldchildlist ) |
---|
| 247 | !--------------------------------------------------------------------------------------------------- |
---|
| 248 | type(Agrif_Grid), pointer, intent(inout) :: g !< Pointer on the current grid |
---|
| 249 | type(Agrif_Grid_List), intent(in) :: oldchildlist |
---|
| 250 | ! |
---|
| 251 | type(Agrif_PGrid), pointer :: parcours ! Pointer for the recursive procedure |
---|
| 252 | real :: g_eps, eps, oldgrid_eps |
---|
| 253 | integer :: out |
---|
| 254 | integer :: iii |
---|
| 255 | ! |
---|
| 256 | out = 0 |
---|
| 257 | ! |
---|
| 258 | parcours => oldchildlist % first |
---|
| 259 | ! |
---|
| 260 | do while (associated(parcours)) |
---|
| 261 | ! |
---|
| 262 | if ((.NOT. g % fixed) .AND. (parcours % gr %oldgrid)) then |
---|
| 263 | ! |
---|
| 264 | g_eps = huge(1.) |
---|
| 265 | oldgrid_eps = huge(1.) |
---|
| 266 | do iii = 1,Agrif_Probdim |
---|
| 267 | g_eps = min(g_eps,g % Agrif_dx(iii)) |
---|
| 268 | oldgrid_eps = min(oldgrid_eps, parcours % gr % Agrif_dx(iii)) |
---|
| 269 | enddo |
---|
| 270 | ! |
---|
| 271 | eps = min(g_eps,oldgrid_eps)/100. |
---|
| 272 | ! |
---|
| 273 | do iii = 1,Agrif_Probdim |
---|
| 274 | if (g % Agrif_dx(iii) < (parcours % gr % Agrif_dx(iii) - eps)) then |
---|
| 275 | parcours => parcours % next |
---|
| 276 | out = 1 |
---|
| 277 | exit |
---|
| 278 | endif |
---|
| 279 | enddo |
---|
| 280 | ! |
---|
| 281 | if ( out == 1 ) cycle |
---|
| 282 | ! |
---|
| 283 | call Agrif_CopyFromOld(g,parcours%gr) |
---|
| 284 | ! |
---|
| 285 | endif |
---|
| 286 | ! |
---|
| 287 | call Agrif_CopyFromOld_All(g,parcours % gr % child_list) |
---|
| 288 | ! |
---|
| 289 | parcours => parcours % next |
---|
| 290 | ! |
---|
| 291 | enddo |
---|
| 292 | !--------------------------------------------------------------------------------------------------- |
---|
| 293 | end subroutine Agrif_CopyFromOld_All |
---|
| 294 | !=================================================================================================== |
---|
| 295 | ! |
---|
| 296 | !=================================================================================================== |
---|
| 297 | ! subroutine Agrif_CopyFromOld |
---|
| 298 | ! |
---|
| 299 | !> Call to the Agrif_Copy procedure. |
---|
| 300 | !--------------------------------------------------------------------------------------------------- |
---|
| 301 | subroutine Agrif_CopyFromOld ( new_gr, old_gr ) |
---|
| 302 | !--------------------------------------------------------------------------------------------------- |
---|
| 303 | type(Agrif_Grid), pointer, intent(inout) :: new_gr !< Pointer on the current grid |
---|
| 304 | type(Agrif_Grid), pointer, intent(inout) :: old_gr !< Pointer on an old grid |
---|
| 305 | ! |
---|
| 306 | type(Agrif_Variable), pointer :: new_var |
---|
| 307 | type(Agrif_Variable), pointer :: old_var |
---|
| 308 | integer :: i |
---|
| 309 | ! |
---|
| 310 | do i = 1,Agrif_NbVariables(0) |
---|
| 311 | if ( Agrif_Mygrid % tabvars(i) % restore ) then |
---|
| 312 | old_var => old_gr % tabvars(i) |
---|
| 313 | new_var => new_gr % tabvars(i) |
---|
| 314 | call Agrif_Copy(new_gr, old_gr, new_var, old_var ) |
---|
| 315 | endif |
---|
| 316 | enddo |
---|
| 317 | !--------------------------------------------------------------------------------------------------- |
---|
| 318 | end subroutine Agrif_CopyFromOld |
---|
| 319 | !=================================================================================================== |
---|
| 320 | ! |
---|
| 321 | !=================================================================================================== |
---|
| 322 | ! subroutine Agrif_CopyFromOld_AllOneVar |
---|
| 323 | ! |
---|
| 324 | !> Called in Agrif_Clustering#Agrif_Init_Hierarchy. |
---|
| 325 | !--------------------------------------------------------------------------------------------------- |
---|
| 326 | recursive subroutine Agrif_CopyFromOld_AllOneVar ( g, oldchildlist, indic ) |
---|
| 327 | !--------------------------------------------------------------------------------------------------- |
---|
| 328 | type(Agrif_Grid), pointer, intent(inout) :: g !< Pointer on the current grid |
---|
| 329 | type(Agrif_Grid_List), intent(in) :: oldchildlist |
---|
| 330 | integer, intent(in) :: indic |
---|
| 331 | ! |
---|
| 332 | type(Agrif_PGrid), pointer :: parcours ! Pointer for the recursive procedure |
---|
| 333 | real :: g_eps,eps,oldgrid_eps |
---|
| 334 | integer :: out |
---|
| 335 | integer :: iii |
---|
| 336 | ! |
---|
| 337 | out = 0 |
---|
| 338 | ! |
---|
| 339 | parcours => oldchildlist % first |
---|
| 340 | ! |
---|
| 341 | do while (associated(parcours)) |
---|
| 342 | ! |
---|
| 343 | if ((.NOT. g % fixed) .AND. (parcours % gr %oldgrid)) then |
---|
| 344 | ! |
---|
| 345 | g_eps = huge(1.) |
---|
| 346 | oldgrid_eps = huge(1.) |
---|
| 347 | do iii = 1,Agrif_Probdim |
---|
| 348 | g_eps = min(g_eps,g % Agrif_dx(iii)) |
---|
| 349 | oldgrid_eps = min(oldgrid_eps,parcours % gr % Agrif_dx(iii)) |
---|
| 350 | enddo |
---|
| 351 | ! |
---|
| 352 | eps = min(g_eps,oldgrid_eps)/100. |
---|
| 353 | ! |
---|
| 354 | do iii = 1,Agrif_Probdim |
---|
| 355 | if (g % Agrif_dx(iii) < (parcours % gr % Agrif_dx(iii) - eps)) then |
---|
| 356 | parcours => parcours % next |
---|
| 357 | out = 1 |
---|
| 358 | exit |
---|
| 359 | endif |
---|
| 360 | enddo |
---|
| 361 | |
---|
| 362 | if ( out == 1 ) cycle |
---|
| 363 | ! |
---|
| 364 | call Agrif_CopyFromOldOneVar(g,parcours%gr,indic) |
---|
| 365 | ! |
---|
| 366 | endif |
---|
| 367 | ! |
---|
| 368 | call Agrif_CopyFromOld_AllOneVar(g,parcours%gr % child_list,indic) |
---|
| 369 | ! |
---|
| 370 | parcours => parcours % next |
---|
| 371 | ! |
---|
| 372 | enddo |
---|
| 373 | !--------------------------------------------------------------------------------------------------- |
---|
| 374 | end subroutine Agrif_CopyFromOld_AllOneVar |
---|
| 375 | !=================================================================================================== |
---|
| 376 | ! |
---|
| 377 | !=================================================================================================== |
---|
| 378 | ! subroutine Agrif_CopyFromOldOneVar |
---|
| 379 | ! |
---|
| 380 | !> Call to Agrif_Copy |
---|
| 381 | !--------------------------------------------------------------------------------------------------- |
---|
| 382 | subroutine Agrif_CopyFromOldOneVar ( new_gr, old_gr, indic ) |
---|
| 383 | !--------------------------------------------------------------------------------------------------- |
---|
| 384 | type(Agrif_Grid), pointer, intent(inout) :: new_gr !< Pointer on the current grid |
---|
| 385 | type(Agrif_Grid), pointer, intent(in) :: old_gr !< Pointer on an old grid |
---|
| 386 | integer, intent(in) :: indic |
---|
| 387 | ! |
---|
| 388 | type(Agrif_Variable), pointer :: new_var, old_var |
---|
| 389 | ! |
---|
| 390 | new_var => Agrif_Search_Variable(new_gr, -indic) |
---|
| 391 | old_var => Agrif_Search_Variable(old_gr, -indic) |
---|
| 392 | ! |
---|
| 393 | call Agrif_array_allocate(new_var,new_var%lb,new_var%ub,new_var%nbdim) |
---|
| 394 | call Agrif_Copy(new_gr, old_gr, new_var,old_var) |
---|
| 395 | !--------------------------------------------------------------------------------------------------- |
---|
| 396 | end subroutine Agrif_CopyFromOldOneVar |
---|
| 397 | !=================================================================================================== |
---|
| 398 | ! |
---|
| 399 | !=================================================================================================== |
---|
| 400 | ! subroutine Agrif_Copy |
---|
| 401 | !--------------------------------------------------------------------------------------------------- |
---|
| 402 | subroutine Agrif_Copy ( new_gr, old_gr, new_var, old_var ) |
---|
| 403 | !--------------------------------------------------------------------------------------------------- |
---|
| 404 | type(Agrif_Grid), pointer, intent(in) :: new_gr !< Pointer on the current grid |
---|
| 405 | type(Agrif_Grid), pointer, intent(in) :: old_gr !< Pointer on an old grid |
---|
| 406 | type(Agrif_Variable), pointer, intent(inout) :: new_var |
---|
| 407 | type(Agrif_Variable), pointer, intent(in) :: old_var |
---|
| 408 | ! |
---|
| 409 | type(Agrif_Variable), pointer :: root ! Pointer on the variable of the root grid |
---|
| 410 | integer :: nbdim ! Number of dimensions of the current grid |
---|
| 411 | integer, dimension(6) :: pttabnew ! Indexes of the first point in the domain |
---|
| 412 | integer, dimension(6) :: petabnew ! Indexes of the first point in the domain |
---|
| 413 | integer, dimension(6) :: pttabold ! Indexes of the first point in the domain |
---|
| 414 | integer, dimension(6) :: petabold ! Indexes of the first point in the domain |
---|
| 415 | integer, dimension(6) :: nbtabold ! Number of cells in each direction |
---|
| 416 | integer, dimension(6) :: nbtabnew ! Number of cells in each direction |
---|
| 417 | real, dimension(6) :: snew,sold |
---|
| 418 | real, dimension(6) :: dsnew,dsold |
---|
| 419 | real :: eps |
---|
| 420 | integer :: n |
---|
| 421 | ! |
---|
| 422 | root => new_var % root_var |
---|
| 423 | nbdim = root % nbdim |
---|
| 424 | ! |
---|
| 425 | do n = 1,nbdim |
---|
| 426 | ! |
---|
| 427 | select case(root % interptab(n)) |
---|
| 428 | ! |
---|
| 429 | case('x') |
---|
| 430 | ! |
---|
| 431 | pttabnew(n) = root % point(n) |
---|
| 432 | pttabold(n) = root % point(n) |
---|
| 433 | snew(n) = new_gr % Agrif_x(1) |
---|
| 434 | sold(n) = old_gr % Agrif_x(1) |
---|
| 435 | dsnew(n) = new_gr % Agrif_dx(1) |
---|
| 436 | dsold(n) = old_gr % Agrif_dx(1) |
---|
| 437 | ! |
---|
| 438 | if (root % posvar(n) == 1) then |
---|
| 439 | petabnew(n) = pttabnew(n) + new_gr % nb(1) |
---|
| 440 | petabold(n) = pttabold(n) + old_gr % nb(1) |
---|
| 441 | else |
---|
| 442 | petabnew(n) = pttabnew(n) + new_gr % nb(1) - 1 |
---|
| 443 | petabold(n) = pttabold(n) + old_gr % nb(1) - 1 |
---|
| 444 | snew(n) = snew(n) + dsnew(n)/2. |
---|
| 445 | sold(n) = sold(n) + dsold(n)/2. |
---|
| 446 | endif |
---|
| 447 | ! |
---|
| 448 | case('y') |
---|
| 449 | ! |
---|
| 450 | pttabnew(n) = root % point(n) |
---|
| 451 | pttabold(n) = root % point(n) |
---|
| 452 | snew(n) = new_gr % Agrif_x(2) |
---|
| 453 | sold(n) = old_gr % Agrif_x(2) |
---|
| 454 | dsnew(n) = new_gr % Agrif_dx(2) |
---|
| 455 | dsold(n) = old_gr % Agrif_dx(2) |
---|
| 456 | ! |
---|
| 457 | if (root % posvar(n) == 1) then |
---|
| 458 | petabnew(n) = pttabnew(n) + new_gr % nb(2) |
---|
| 459 | petabold(n) = pttabold(n) + old_gr % nb(2) |
---|
| 460 | else |
---|
| 461 | petabnew(n) = pttabnew(n) + new_gr % nb(2) - 1 |
---|
| 462 | petabold(n) = pttabold(n) + old_gr % nb(2) - 1 |
---|
| 463 | snew(n) = snew(n) + dsnew(n)/2. |
---|
| 464 | sold(n) = sold(n) + dsold(n)/2. |
---|
| 465 | endif |
---|
| 466 | ! |
---|
| 467 | case('z') |
---|
| 468 | ! |
---|
| 469 | pttabnew(n) = root % point(n) |
---|
| 470 | pttabold(n) = root % point(n) |
---|
| 471 | snew(n) = new_gr % Agrif_x(3) |
---|
| 472 | sold(n) = old_gr % Agrif_x(3) |
---|
| 473 | dsnew(n) = new_gr % Agrif_dx(3) |
---|
| 474 | dsold(n) = old_gr % Agrif_dx(3) |
---|
| 475 | ! |
---|
| 476 | if (root % posvar(n) == 1) then |
---|
| 477 | petabnew(n) = pttabnew(n) + new_gr % nb(3) |
---|
| 478 | petabold(n) = pttabold(n) + old_gr % nb(3) |
---|
| 479 | else |
---|
| 480 | petabnew(n) = pttabnew(n) + new_gr % nb(3) - 1 |
---|
| 481 | petabold(n) = pttabold(n) + old_gr % nb(3) - 1 |
---|
| 482 | snew(n) = snew(n) + dsnew(n)/2. |
---|
| 483 | sold(n) = sold(n) + dsold(n)/2. |
---|
| 484 | endif |
---|
| 485 | ! |
---|
| 486 | case('N') |
---|
| 487 | ! |
---|
| 488 | call Agrif_get_var_bounds(new_var,pttabnew(n),petabnew(n),n) |
---|
| 489 | ! |
---|
| 490 | pttabold(n) = pttabnew(n) |
---|
| 491 | petabold(n) = petabnew(n) |
---|
| 492 | snew(n) = 0. |
---|
| 493 | sold(n) = 0. |
---|
| 494 | dsnew(n) = 1. |
---|
| 495 | dsold(n) = 1. |
---|
| 496 | ! |
---|
| 497 | end select |
---|
| 498 | ! |
---|
| 499 | enddo |
---|
| 500 | ! |
---|
| 501 | do n = 1,nbdim |
---|
| 502 | nbtabnew(n) = petabnew(n) - pttabnew(n) |
---|
| 503 | nbtabold(n) = petabold(n) - pttabold(n) |
---|
| 504 | enddo |
---|
| 505 | ! |
---|
| 506 | eps = min(minval(dsnew(1:nbdim)),minval(dsold(1:nbdim))) / 100. |
---|
| 507 | ! |
---|
| 508 | do n = 1,nbdim |
---|
| 509 | if (snew(n) > (sold(n)+dsold(n)*nbtabold(n)+eps)) return |
---|
| 510 | if ((snew(n)+dsnew(n)*nbtabnew(n)-eps) < sold(n)) return |
---|
| 511 | enddo |
---|
| 512 | ! |
---|
| 513 | call Agrif_CopynD(new_var,old_var,pttabold,petabold,pttabnew,petabnew, & |
---|
| 514 | sold,snew,dsold,dsnew,nbdim) |
---|
| 515 | !--------------------------------------------------------------------------------------------------- |
---|
| 516 | end subroutine Agrif_Copy |
---|
| 517 | !=================================================================================================== |
---|
| 518 | ! |
---|
| 519 | !=================================================================================================== |
---|
| 520 | ! subroutine Agrif_CopynD |
---|
| 521 | ! |
---|
| 522 | !> Copy from the nD variable Old_Var the nD variable New_Var |
---|
| 523 | !--------------------------------------------------------------------------------------------------- |
---|
| 524 | subroutine Agrif_CopynD ( new_var, old_var, pttabold, petabold, pttabnew, petabnew, & |
---|
| 525 | sold, snew, dsold, dsnew, nbdim ) |
---|
| 526 | !--------------------------------------------------------------------------------------------------- |
---|
| 527 | type(Agrif_Variable), pointer, intent(inout) :: new_var |
---|
| 528 | type(Agrif_Variable), pointer, intent(in) :: old_var |
---|
| 529 | integer, dimension(nbdim), intent(in) :: pttabnew |
---|
| 530 | integer, dimension(nbdim), intent(in) :: petabnew |
---|
| 531 | integer, dimension(nbdim), intent(in) :: pttabold |
---|
| 532 | integer, dimension(nbdim), intent(in) :: petabold |
---|
| 533 | real, dimension(nbdim), intent(in) :: snew, sold |
---|
| 534 | real, dimension(nbdim), intent(in) :: dsnew,dsold |
---|
| 535 | integer, intent(in) :: nbdim |
---|
| 536 | ! |
---|
| 537 | integer :: i,j,k,l,m,n,i0,j0,k0,l0,m0,n0 |
---|
| 538 | ! |
---|
| 539 | real, dimension(nbdim) :: dim_gmin, dim_gmax |
---|
| 540 | real, dimension(nbdim) :: dim_newmin, dim_newmax |
---|
| 541 | real, dimension(nbdim) :: dim_min |
---|
| 542 | integer, dimension(nbdim) :: ind_gmin,ind_newmin, ind_newmax |
---|
| 543 | ! |
---|
| 544 | do i = 1,nbdim |
---|
| 545 | ! |
---|
| 546 | dim_gmin(i) = sold(i) |
---|
| 547 | dim_gmax(i) = sold(i) + (petabold(i)-pttabold(i)) * dsold(i) |
---|
| 548 | ! |
---|
| 549 | dim_newmin(i) = snew(i) |
---|
| 550 | dim_newmax(i) = snew(i) + (petabnew(i)-pttabnew(i)) * dsnew(i) |
---|
| 551 | ! |
---|
| 552 | enddo |
---|
| 553 | ! |
---|
| 554 | do i = 1,nbdim |
---|
| 555 | if (dim_gmax(i) < dim_newmin(i)) return |
---|
| 556 | if (dim_gmin(i) > dim_newmax(i)) return |
---|
| 557 | enddo |
---|
| 558 | ! |
---|
| 559 | do i = 1,nbdim |
---|
| 560 | ! |
---|
| 561 | ind_newmin(i) = pttabnew(i) - floor(-(max(dim_gmin(i),dim_newmin(i))-dim_newmin(i))/dsnew(i)) |
---|
| 562 | dim_min(i) = snew(i) + (ind_newmin(i)-pttabnew(i))*dsnew(i) |
---|
| 563 | ind_gmin(i) = pttabold(i) + nint((dim_min(i)-dim_gmin(i))/dsold(i)) |
---|
| 564 | ind_newmax(i) = pttabnew(i) + int((min(dim_gmax(i),dim_newmax(i))-dim_newmin(i))/dsnew(i)) |
---|
| 565 | ! |
---|
| 566 | enddo |
---|
| 567 | ! |
---|
| 568 | select case (nbdim) |
---|
| 569 | ! |
---|
| 570 | case (1) |
---|
| 571 | i0 = ind_gmin(1) |
---|
| 572 | do i = ind_newmin(1),ind_newmax(1) |
---|
| 573 | new_var % array1(i) = old_var % array1(i0) |
---|
| 574 | new_var % restore1D(i) = 1 |
---|
| 575 | i0 = i0 + int(dsnew(1)/dsold(1)) |
---|
| 576 | enddo |
---|
| 577 | ! |
---|
| 578 | case (2) |
---|
| 579 | i0 = ind_gmin(1) ; do i = ind_newmin(1),ind_newmax(1) |
---|
| 580 | j0 = ind_gmin(2) ; do j = ind_newmin(2),ind_newmax(2) |
---|
| 581 | new_var % array2(i,j) = old_var % array2(i0,j0) |
---|
| 582 | new_var % restore2D(i,j) = 1 |
---|
| 583 | j0 = j0 + int(dsnew(2)/dsold(2)) |
---|
| 584 | enddo |
---|
| 585 | i0 = i0 + int(dsnew(1)/dsold(1)) |
---|
| 586 | enddo |
---|
| 587 | ! |
---|
| 588 | case (3) |
---|
| 589 | i0 = ind_gmin(1) ; do i = ind_newmin(1),ind_newmax(1) |
---|
| 590 | j0 = ind_gmin(2) ; do j = ind_newmin(2),ind_newmax(2) |
---|
| 591 | k0 = ind_gmin(3) ; do k = ind_newmin(3),ind_newmax(3) |
---|
| 592 | new_var % array3(i,j,k) = old_var % array3(i0,j0,k0) |
---|
| 593 | new_var % restore3D(i,j,k) = 1 |
---|
| 594 | k0 = k0 + int(dsnew(3)/dsold(3)) |
---|
| 595 | enddo |
---|
| 596 | j0 = j0 + int(dsnew(2)/dsold(2)) |
---|
| 597 | enddo |
---|
| 598 | i0 = i0 + int(dsnew(1)/dsold(1)) |
---|
| 599 | enddo |
---|
| 600 | ! |
---|
| 601 | case (4) |
---|
| 602 | i0 = ind_gmin(1) ; do i = ind_newmin(1),ind_newmax(1) |
---|
| 603 | j0 = ind_gmin(2) ; do j = ind_newmin(2),ind_newmax(2) |
---|
| 604 | k0 = ind_gmin(3) ; do k = ind_newmin(3),ind_newmax(3) |
---|
| 605 | l0 = ind_gmin(4) ; do l = ind_newmin(4),ind_newmax(4) |
---|
| 606 | new_var % array4(i,j,k,l) = old_var % array4(i0,j0,k0,l0) |
---|
| 607 | new_var % restore4D(i,j,k,l) = 1 |
---|
| 608 | l0 = l0 + int(dsnew(4)/dsold(4)) |
---|
| 609 | enddo |
---|
| 610 | k0 = k0 + int(dsnew(3)/dsold(3)) |
---|
| 611 | enddo |
---|
| 612 | j0 = j0 + int(dsnew(2)/dsold(2)) |
---|
| 613 | enddo |
---|
| 614 | i0 = i0 + int(dsnew(1)/dsold(1)) |
---|
| 615 | enddo |
---|
| 616 | ! |
---|
| 617 | case (5) |
---|
| 618 | i0 = ind_gmin(1) ; do i = ind_newmin(1),ind_newmax(1) |
---|
| 619 | j0 = ind_gmin(2) ; do j = ind_newmin(2),ind_newmax(2) |
---|
| 620 | k0 = ind_gmin(3) ; do k = ind_newmin(3),ind_newmax(3) |
---|
| 621 | l0 = ind_gmin(4) ; do l = ind_newmin(4),ind_newmax(4) |
---|
| 622 | m0 = ind_gmin(5) ; do m = ind_newmin(5),ind_newmax(5) |
---|
| 623 | new_var % array5(i,j,k,l,m) = old_var % array5(i0,j0,k0,l0,m0) |
---|
| 624 | new_var % restore5D(i,j,k,l,m) = 1 |
---|
| 625 | m0 = m0 + int(dsnew(5)/dsold(5)) |
---|
| 626 | enddo |
---|
| 627 | l0 = l0 + int(dsnew(4)/dsold(4)) |
---|
| 628 | enddo |
---|
| 629 | k0 = k0 + int(dsnew(3)/dsold(3)) |
---|
| 630 | enddo |
---|
| 631 | j0 = j0 + int(dsnew(2)/dsold(2)) |
---|
| 632 | enddo |
---|
| 633 | i0 = i0 + int(dsnew(1)/dsold(1)) |
---|
| 634 | enddo |
---|
| 635 | ! |
---|
| 636 | case (6) |
---|
| 637 | i0 = ind_gmin(1) ; do i = ind_newmin(1),ind_newmax(1) |
---|
| 638 | j0 = ind_gmin(2) ; do j = ind_newmin(2),ind_newmax(2) |
---|
| 639 | k0 = ind_gmin(3) ; do k = ind_newmin(3),ind_newmax(3) |
---|
| 640 | l0 = ind_gmin(4) ; do l = ind_newmin(4),ind_newmax(4) |
---|
| 641 | m0 = ind_gmin(5) ; do m = ind_newmin(5),ind_newmax(5) |
---|
| 642 | n0 = ind_gmin(6) ; do n = ind_newmin(6),ind_newmax(6) |
---|
| 643 | new_var % array6(i,j,k,l,m,n) = old_var % array6(i0,j0,k0,l0,m0,n0) |
---|
| 644 | new_var % restore6D(i,j,k,l,m,n) = 1 |
---|
| 645 | n0 = n0 + int(dsnew(6)/dsold(6)) |
---|
| 646 | enddo |
---|
| 647 | m0 = m0 + int(dsnew(5)/dsold(5)) |
---|
| 648 | enddo |
---|
| 649 | l0 = l0 + int(dsnew(4)/dsold(4)) |
---|
| 650 | enddo |
---|
| 651 | k0 = k0 + int(dsnew(3)/dsold(3)) |
---|
| 652 | enddo |
---|
| 653 | j0 = j0 + int(dsnew(2)/dsold(2)) |
---|
| 654 | enddo |
---|
| 655 | i0 = i0 + int(dsnew(1)/dsold(1)) |
---|
| 656 | enddo |
---|
| 657 | ! |
---|
| 658 | end select |
---|
| 659 | !--------------------------------------------------------------------------------------------------- |
---|
| 660 | end subroutine Agrif_CopynD |
---|
| 661 | !=================================================================================================== |
---|
| 662 | ! |
---|
| 663 | end module Agrif_Save |
---|