Rocstar  1.0
Rocstar multiphysics simulation application
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
RFLU_ModSymmetryPeriodic.F90
Go to the documentation of this file.
1 ! *********************************************************************
2 ! * Rocstar Simulation Suite *
3 ! * Copyright@2015, Illinois Rocstar LLC. All rights reserved. *
4 ! * *
5 ! * Illinois Rocstar LLC *
6 ! * Champaign, IL *
7 ! * www.illinoisrocstar.com *
8 ! * sales@illinoisrocstar.com *
9 ! * *
10 ! * License: See LICENSE file in top level of distribution package or *
11 ! * http://opensource.org/licenses/NCSA *
12 ! *********************************************************************
13 ! *********************************************************************
14 ! * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, *
15 ! * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES *
16 ! * OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND *
17 ! * NONINFRINGEMENT. IN NO EVENT SHALL THE CONTRIBUTORS OR *
18 ! * COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER *
19 ! * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, *
20 ! * Arising FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE *
21 ! * USE OR OTHER DEALINGS WITH THE SOFTWARE. *
22 ! *********************************************************************
23 ! ******************************************************************************
24 !
25 ! Purpose: Collection of routines for symmetry and periodic patch operations.
26 !
27 ! Description: None.
28 !
29 ! Notes: None.
30 !
31 ! ******************************************************************************
32 !
33 ! $Id: RFLU_ModSymmetryPeriodic.F90,v 1.8 2008/12/06 08:44:24 mtcampbe Exp $
34 !
35 ! Copyright: (c) 2006 by the University of Illinois
36 !
37 ! ******************************************************************************
38 
40 
41  USE modparameters
42  USE moddatatypes
43  USE modglobal, ONLY: t_global
44  USE modgrid, ONLY: t_grid
45  USE modbndpatch, ONLY: t_patch
46  USE moddatastruct, ONLY: t_region
47  USE modborder, ONLY: t_border
48  USE moderror
49  USE modmpi
50 
51  USE modsortsearch
52 
53  IMPLICIT NONE
54 
55 ! ******************************************************************************
56 ! Definitions and declarations
57 ! ******************************************************************************
58 
59 ! ==============================================================================
60 ! Private data
61 ! ==============================================================================
62 
63  CHARACTER(CHRLEN), PRIVATE :: &
64  RCSIdentString = '$RCSfile: RFLU_ModSymmetryPeriodic.F90,v $ $Revision: 1.8 $'
65 
66 ! ==============================================================================
67 ! Public functions
68 ! ==============================================================================
69 
70  PUBLIC :: rflu_sype_addvirtualcells, &
83 
84 ! ==============================================================================
85 ! Private functions
86 ! ==============================================================================
87 
89 
90 
91 ! ******************************************************************************
92 ! Routines
93 ! ******************************************************************************
94 
95  CONTAINS
96 
97 
98 
99 
100 
101 
102 ! ******************************************************************************
103 !
104 ! Purpose: Add virtual cells arising from sype patches.
105 !
106 ! Description: None.
107 !
108 ! Input:
109 ! pRegion Pointer to region
110 !
111 ! Output: None.
112 !
113 ! Notes: None.
114 !
115 ! ******************************************************************************
116 
117 SUBROUTINE rflu_sype_addvirtualcells(pRegion)
118 
120 
121  IMPLICIT NONE
122 
123 ! ******************************************************************************
124 ! Declarations and definitions
125 ! ******************************************************************************
126 
127 ! ==============================================================================
128 ! Arguments
129 ! ==============================================================================
130 
131  TYPE(t_region), POINTER :: pregion
132 
133 ! ==============================================================================
134 ! Locals
135 ! ==============================================================================
136 
137  INTEGER :: errorflag,i,icg,icg2,icl,ict,ilayer,iloc,ipatch,ireg,j,key, &
138  ncellsvirt,ncellsvirtmax,nlayers
139  INTEGER, DIMENSION(:), ALLOCATABLE :: vc
140  TYPE(t_global), POINTER :: global
141  TYPE(t_grid), POINTER :: pgrid
142  TYPE(t_patch), POINTER :: ppatch,ppatchrelated
143 
144 ! ******************************************************************************
145 ! Start
146 ! ******************************************************************************
147 
148  global => pregion%global
149 
150  CALL registerfunction(global,'RFLU_SYPE_AddVirtualCells',&
151  'RFLU_ModSymmetryPeriodic.F90')
152 
153  IF ( global%verbLevel > verbose_none ) THEN
154  WRITE(stdout,'(A,1X,A)') solver_name, &
155  'Adding s/p virtual cells...'
156  END IF ! global%verbLevel
157 
158  IF ( global%verbLevel > verbose_low ) THEN
159  WRITE(stdout,'(A,3X,A,1X,I5.5)') solver_name,'Global region:', &
160  pregion%iRegionGlobal
161  END IF ! global%verbLevel
162 
163 ! ******************************************************************************
164 ! Set pointers
165 ! ******************************************************************************
166 
167  pgrid => pregion%grid
168 
169 ! ******************************************************************************
170 ! Allocate temporary for virtual cells
171 ! ******************************************************************************
172 
173  ncellsvirtmax = pgrid%nCellsMax - pgrid%nCells
174 
175  ALLOCATE(vc(ncellsvirtmax),stat=errorflag)
176  global%error = errorflag
177  IF ( global%error /= err_none ) THEN
178  CALL errorstop(global,err_allocate,__line__,'vc')
179  END IF ! global%error
180 
181 ! ******************************************************************************
182 ! Add virtual cells
183 ! ******************************************************************************
184 
185  DO ipatch = 1,pgrid%nPatches
186  ppatch => pregion%patches(ipatch)
187 
188  SELECT CASE ( ppatch%bcType )
189  CASE ( bc_symmetry )
190  IF ( global%verbLevel > verbose_low ) THEN
191  WRITE(stdout,'(A,5X,A,2(1X,I3))') solver_name,'Patch:',ipatch, &
192  ppatch%iPatchGlobal
193  END IF ! global%verbLevel
194 
195  SELECT CASE ( pregion%mixtInput%spaceOrder )
196  CASE ( 1 )
197  CALL rflu_sype_addvirtualcellsinv1(pregion,ppatch,vc, &
198  ncellsvirtmax,ncellsvirt)
199  CASE ( 2 )
200  CALL rflu_sype_addvirtualcellsinv2(pregion,ppatch,vc, &
201  ncellsvirtmax,ncellsvirt)
202  CASE default
203  CALL errorstop(global,err_reached_default,__line__)
204  END SELECT ! pRegion%mixtInput%spaceOrder
205 
206  CALL rflu_sype_augmentcelllists(pregion,ppatch,vc,ncellsvirt)
207  CASE ( bc_periodic )
208  ppatchrelated => pregion%patches(ppatch%iPatchRelated)
209 
210  IF ( global%verbLevel > verbose_low ) THEN
211  WRITE(stdout,'(A,5X,A,2(1X,I3))') solver_name,'Patch:',ipatch, &
212  ppatch%iPatchGlobal
213  END IF ! global%verbLevel
214 
215  SELECT CASE ( pregion%mixtInput%spaceOrder )
216  CASE ( 1 )
217  CALL rflu_sype_addvirtualcellsinv1(pregion,ppatchrelated,vc, &
218  ncellsvirtmax,ncellsvirt)
219  CASE ( 2 )
220  CALL rflu_sype_addvirtualcellsinv2(pregion,ppatchrelated,vc, &
221  ncellsvirtmax,ncellsvirt)
222  CASE default
223  CALL errorstop(global,err_reached_default,__line__)
224  END SELECT ! pRegion%mixtInput%spaceOrder
225 
226  CALL rflu_sype_augmentcelllists(pregion,ppatch,vc,ncellsvirt)
227  END SELECT ! pPatch%bcType
228  END DO ! iPatch
229 
230 ! ******************************************************************************
231 ! Deallocate temporary memory for virtual cells
232 ! ******************************************************************************
233 
234  DEALLOCATE(vc,stat=errorflag)
235  global%error = errorflag
236  IF ( global%error /= err_none ) THEN
237  CALL errorstop(global,err_deallocate,__line__,'virtCells')
238  END IF ! global%error
239 
240 ! ******************************************************************************
241 ! Write information about numbers of virtual cells
242 ! ******************************************************************************
243 
244  IF ( global%verbLevel > verbose_low ) THEN
245  WRITE(stdout,'(A,3X,A)') solver_name,'Virtual cell statistics:'
246  WRITE(stdout,'(A,5X,A,1X,I6)') solver_name,'Tetrahedra:', &
247  pgrid%nTetsTot-pgrid%nTets
248  WRITE(stdout,'(A,5X,A,1X,I6)') solver_name,'Hexahedra: ', &
249  pgrid%nHexsTot-pgrid%nHexs
250  WRITE(stdout,'(A,5X,A,1X,I6)') solver_name,'Prisms: ', &
251  pgrid%nPrisTot-pgrid%nPris
252  WRITE(stdout,'(A,5X,A,1X,I6)') solver_name,'Pyramids: ', &
253  pgrid%nPyrsTot-pgrid%nPyrs
254  END IF ! global%verbLevel
255 
256 ! ******************************************************************************
257 ! End
258 ! ******************************************************************************
259 
260  IF ( global%verbLevel > verbose_none ) THEN
261  WRITE(stdout,'(A,1X,A)') solver_name, &
262  'Adding s/p virtual cells done.'
263  END IF ! global%verbLevel
264 
265  CALL deregisterfunction(global)
266 
267 END SUBROUTINE rflu_sype_addvirtualcells
268 
269 
270 
271 
272 
273 
274 
275 ! ******************************************************************************
276 !
277 ! Purpose: Add virtual cells for inviscid fluxes based on first-order scheme.
278 !
279 ! Description: None.
280 !
281 ! Input:
282 ! pRegion Pointer to region
283 ! pRegionSerial Pointer to serial region
284 ! vc List of virtual cells (empty)
285 ! nCellsVirtMax Maximum allowable number of virtual cells
286 !
287 ! Output:
288 ! vc List of virtual cells
289 ! nCellsVirt Number of virtual cells
290 !
291 ! Notes: None.
292 !
293 ! ******************************************************************************
294 
295 SUBROUTINE rflu_sype_addvirtualcellsinv1(pRegion,pPatch,vc,nCellsVirtMax, &
296  ncellsvirt)
297 
298  IMPLICIT NONE
299 
300 ! ******************************************************************************
301 ! Declarations and definitions
302 ! ******************************************************************************
303 
304 ! ==============================================================================
305 ! Arguments
306 ! ==============================================================================
307 
308  INTEGER, INTENT(IN) :: ncellsvirtmax
309  INTEGER, INTENT(OUT) :: ncellsvirt
310  INTEGER, INTENT(INOUT) :: vc(ncellsvirtmax)
311  TYPE(t_patch), POINTER :: ppatch
312  TYPE(t_region), POINTER :: pregion
313 
314 ! ==============================================================================
315 ! Locals
316 ! ==============================================================================
317 
318  INTEGER :: errorflag,ifl
319  TYPE(t_global), POINTER :: global
320 
321 ! ******************************************************************************
322 ! Start
323 ! ******************************************************************************
324 
325  global => pregion%global
326 
327  CALL registerfunction(global,'RFLU_SYPE_AddVirtualCellsInv1',&
328  'RFLU_ModSymmetryPeriodic.F90')
329 
330  IF ( global%myProcid == masterproc .AND. &
331  global%verbLevel > verbose_none ) THEN
332  WRITE(stdout,'(A,1X,A)') solver_name, &
333  'Adding s/p virtual cells for inviscid first-order stencil...'
334  END IF ! global%myProcid
335 
336 ! ******************************************************************************
337 ! Set pointers and initialize variables
338 ! ******************************************************************************
339 
340  ncellsvirt = 0
341 
342 ! ******************************************************************************
343 ! Loop over faces on this patch and add adjacent cells
344 ! ******************************************************************************
345 
346  DO ifl = 1,ppatch%nBFaces
347  IF ( ncellsvirt < ncellsvirtmax ) THEN
348  ncellsvirt = ncellsvirt + 1
349  vc(ncellsvirt) = ppatch%bf2c(ifl)
350  ELSE
351  CALL errorstop(global,err_exceed_dimens,__line__,'vc')
352  END IF ! nCellsVirt
353  END DO ! ifl
354 
355 ! ******************************************************************************
356 ! End
357 ! ******************************************************************************
358 
359  IF ( global%myProcid == masterproc .AND. &
360  global%verbLevel > verbose_none ) THEN
361  WRITE(stdout,'(A,1X,A)') solver_name, &
362  'Adding s/p virtual cells for first-order stencil done.'
363  END IF ! global%myProcid
364 
365  CALL deregisterfunction(global)
366 
367 END SUBROUTINE rflu_sype_addvirtualcellsinv1
368 
369 
370 
371 
372 
373 
374 
375 
376 ! ******************************************************************************
377 !
378 ! Purpose: Add virtual cells for inviscid fluxes based on higher-order scheme.
379 !
380 ! Description: Building list of virtual cells proceeds in several steps:
381 ! 1. Build list of vertices from list of actual-virtual faces
382 ! 2. Build list of source cells adjacent to actual-virtual faces which are
383 ! in same region as that for which virtual cells are constructed. One
384 ! layer of source cells is sufficient. If more layers of source cells
385 ! are specified here, this should not change the number of virtual cells.
386 ! 3. Build stencils for source cells and add stencil members to list of
387 ! virtual cells if not in same region
388 ! 4. Loop over existing virtual cells and add layers.
389 !
390 ! Input:
391 ! pRegion Pointer to region
392 ! pRegionSerial Pointer to serial region
393 ! vc List of virtual cells (empty)
394 ! nCellsVirtMax Maximum allowable number of virtual cells
395 !
396 ! Output:
397 ! vc List of virtual cells
398 ! nCellsVirt Number of virtual cells
399 !
400 ! Notes: None.
401 !
402 ! ******************************************************************************
403 
404 SUBROUTINE rflu_sype_addvirtualcellsinv2(pRegion,pPatch,vc,nCellsVirtMax, &
405  ncellsvirt)
406 
407  USE modsortsearch
408 
415 
416  IMPLICIT NONE
417 
418 ! ******************************************************************************
419 ! Declarations and definitions
420 ! ******************************************************************************
421 
422 ! ==============================================================================
423 ! Arguments
424 ! ==============================================================================
425 
426  INTEGER, INTENT(IN) :: ncellsvirtmax
427  INTEGER, INTENT(OUT) :: ncellsvirt
428  INTEGER, INTENT(INOUT) :: vc(ncellsvirtmax)
429  TYPE(t_patch), POINTER :: ppatch
430  TYPE(t_region), POINTER :: pregion
431 
432 ! ==============================================================================
433 ! Locals
434 ! ==============================================================================
435 
436  INTEGER :: errorflag,i,icg,icg2,iflbeg,iflend,ilayer,iloc,j,nlayers, &
437  nvert,scdim,vcnewdim,vcnewdimmax,vcolddim
438  INTEGER, DIMENSION(:), ALLOCATABLE :: avv,sc,vcnew,vcold
439  TYPE(t_grid), POINTER :: pgrid
440  TYPE(t_global), POINTER :: global
441 
442 ! ******************************************************************************
443 ! Start
444 ! ******************************************************************************
445 
446  global => pregion%global
447 
448  CALL registerfunction(global,'RFLU_SYPE_AddVirtualCellsInv2',&
449  'RFLU_ModSymmetryPeriodic.F90')
450 
451  IF ( global%myProcid == masterproc .AND. &
452  global%verbLevel > verbose_none ) THEN
453  WRITE(stdout,'(A,1X,A)') solver_name, &
454  'Adding s/p virtual cells for inviscid higher-order stencil...'
455  END IF ! global%myProcid
456 
457 ! ******************************************************************************
458 ! Set pointers and initialize variables
459 ! ******************************************************************************
460 
461  pgrid => pregion%grid
462 
463  ncellsvirt = 0
464 
465 ! ******************************************************************************
466 ! Based on list of patch vertices, build list of source cells adjacent to
467 ! actual-virtual faces for construction of virtual cells. One layer of source
468 ! cells is sufficient. If more layers of source cells are specified here, this
469 ! should not change the number of virtual cells.
470 ! ******************************************************************************
471 
472  ALLOCATE(sc(ncellsvirtmax),stat=errorflag)
473  global%error = errorflag
474  IF ( global%error /= err_none ) THEN
475  CALL errorstop(global,err_allocate,__line__,'sc')
476  END IF ! global%error
477 
478  nlayers = 1
479 
480  CALL rflu_buildvertcellnghblist(global,pgrid,ppatch%bv,ppatch%nBVert, &
481  nlayers,pregion%iRegionGlobal,sc, &
482  ncellsvirtmax,scdim)
483 
484 ! ******************************************************************************
485 ! Loop over list of source cells, build stencil for each cell, and add cells
486 ! in stencil as virtual cells if not in same region and not already added
487 ! ******************************************************************************
488 
489  CALL rflu_setinfoc2cstencilwrapper(pregion,pregion%mixtInput%spaceOrder-1)
490  CALL rflu_createc2cstencilwrapper(pregion)
491 
492  DO i = 1,scdim
493  icg = sc(i)
494 
495  CALL rflu_buildc2cstencilwrapper(pregion,icg,constr_none)
496 
497  DO j = 1,pgrid%c2cs(icg)%nCellMembs
498  icg2 = pgrid%c2cs(icg)%cellMembs(j)
499 
500  IF ( ncellsvirt > 0 ) THEN
501  CALL binarysearchinteger(vc(1:ncellsvirt),ncellsvirt,icg2,iloc)
502  ELSE
503  iloc = element_not_found
504  END IF ! nCellsVirt
505 
506  IF ( iloc == element_not_found ) THEN
507  IF ( ncellsvirt < ncellsvirtmax ) THEN
508  ncellsvirt = ncellsvirt + 1
509  vc(ncellsvirt) = icg2
510  ELSE
511  CALL errorstop(global,err_exceed_dimens,__line__,'vc')
512  END IF ! nCellsVirt
513 
514  IF ( ncellsvirt > 1 ) THEN
515  CALL quicksortinteger(vc(1:ncellsvirt),ncellsvirt)
516  END IF ! nCellsVirt
517  END IF ! iLoc
518  END DO ! j
519  END DO ! i
520 
521  DEALLOCATE(sc,stat=errorflag)
522  global%error = errorflag
523  IF ( global%error /= err_none ) THEN
524  CALL errorstop(global,err_deallocate,__line__,'sc')
525  END IF ! global%error
526 
527 ! ******************************************************************************
528 ! Loop over current list of virtual cells and add layers of cells. For the
529 ! current scheme, two additional layers are required.
530 ! ******************************************************************************
531 
532 ! ==============================================================================
533 ! Allocate temporary memory
534 ! ==============================================================================
535 
536  vcolddim = ncellsvirt
537  vcnewdimmax = ncellsvirtmax
538 
539  ALLOCATE(vcold(vcolddim),stat=errorflag)
540  global%error = errorflag
541  IF ( global%error /= err_none ) THEN
542  CALL errorstop(global,err_allocate,__line__,'vcOld')
543  END IF ! global%error
544 
545  DO i = 1,vcolddim
546  vcold(i) = vc(i)
547  END DO ! i
548 
549  ALLOCATE(vcnew(vcnewdimmax),stat=errorflag)
550  global%error = errorflag
551  IF ( global%error /= err_none ) THEN
552  CALL errorstop(global,err_allocate,__line__,'vcNew')
553  END IF ! global%error
554 
555 ! ==============================================================================
556 ! Loop over layers
557 ! ==============================================================================
558 
559  nlayers = 1
560 
561  DO ilayer = 1,nlayers
562  vcnewdim = 0
563 
564  DO i = 1,vcolddim
565  icg = vcold(i)
566 
567  CALL rflu_buildc2cstencilwrapper(pregion,icg,constr_none)
568 
569  DO j = 1,pgrid%c2cs(icg)%nCellMembs
570  icg2 = pgrid%c2cs(icg)%cellMembs(j)
571 
572  CALL binarysearchinteger(vc(1:ncellsvirt),ncellsvirt,icg2,iloc)
573 
574  IF ( iloc == element_not_found ) THEN
575  IF ( ncellsvirt < ncellsvirtmax ) THEN
576  ncellsvirt = ncellsvirt + 1
577  vc(ncellsvirt) = icg2
578  ELSE
579  CALL errorstop(global,err_exceed_dimens,__line__,'vc')
580  END IF ! nCellsVirt
581 
582  IF ( vcnewdim < vcnewdimmax ) THEN
583  vcnewdim = vcnewdim + 1
584  vcnew(vcnewdim) = icg2
585  ELSE
586  CALL errorstop(global,err_exceed_dimens,__line__,'vcNew')
587  END IF ! vcNewDim
588 
589  CALL quicksortinteger(vc(1:ncellsvirt),ncellsvirt)
590  END IF ! iLoc
591  END DO ! j
592  END DO ! i
593 
594  DEALLOCATE(vcold,stat=errorflag)
595  global%error = errorflag
596  IF ( global%error /= err_none ) THEN
597  CALL errorstop(global,err_deallocate,__line__,'vcOld')
598  END IF ! global%error
599 
600  IF ( ilayer /= nlayers ) THEN
601  vcolddim = vcnewdim
602 
603  ALLOCATE(vcold(vcolddim),stat=errorflag)
604  global%error = errorflag
605  IF ( global%error /= err_none ) THEN
606  CALL errorstop(global,err_allocate,__line__,'vcOld')
607  END IF ! global%error
608 
609  DO i = 1,vcolddim
610  vcold(i) = vcnew(i)
611  END DO ! i
612  END IF ! iLayer
613  END DO ! iLayer
614 
615 ! ==============================================================================
616 ! Deallocate temporary memory
617 ! ==============================================================================
618 
619  DEALLOCATE(vcnew,stat=errorflag)
620  global%error = errorflag
621  IF ( global%error /= err_none ) THEN
622  CALL errorstop(global,err_deallocate,__line__,'vcNew')
623  END IF ! global%error
624 
625  CALL rflu_destroyc2cstencilwrapper(pregion)
626 
627 ! ******************************************************************************
628 ! End
629 ! ******************************************************************************
630 
631  IF ( global%myProcid == masterproc .AND. &
632  global%verbLevel > verbose_none ) THEN
633  WRITE(stdout,'(A,1X,A)') solver_name, &
634  'Adding s/p virtual cells for inviscid higher-order stencil done.'
635  END IF ! global%myProcid
636 
637  CALL deregisterfunction(global)
638 
639 END SUBROUTINE rflu_sype_addvirtualcellsinv2
640 
641 
642 
643 
644 
645 
646 
647 ! ******************************************************************************
648 !
649 ! Purpose: Add virtual cells to cell and patch data structures.
650 !
651 ! Description: None.
652 !
653 ! Input:
654 ! pRegion Pointer to region
655 ! pRegionSerial Pointer to serial region
656 ! vc List of virtual cells
657 ! nCellsVirt Number of virtual cells
658 !
659 ! Output: None.
660 !
661 ! Notes: None.
662 !
663 ! ******************************************************************************
664 
665 SUBROUTINE rflu_sype_augmentcelllists(pRegion,pPatch,vc,nCellsVirt)
666 
667  USE modtools, ONLY: swapintegers, &
669 
670  IMPLICIT NONE
671 
672 ! ******************************************************************************
673 ! Declarations and definitions
674 ! ******************************************************************************
675 
676 ! ==============================================================================
677 ! Arguments
678 ! ==============================================================================
679 
680  INTEGER, INTENT(IN) :: ncellsvirt
681  INTEGER, INTENT(IN) :: vc(ncellsvirt)
682  TYPE(t_patch), POINTER :: ppatch
683  TYPE(t_region), POINTER :: pregion
684 
685 ! ==============================================================================
686 ! Locals
687 ! ==============================================================================
688 
689  INTEGER :: errorflag,iborder,icg,icg2,icl,ict,icv,ifl,iflbeg,iflend,ifl2, &
690  ilayer,iloc,iloc2,ipatchrelated,ipatch2,ivg,ivgm,ivl,ivlm,ivl2, &
691  j,nlayers,nvert,nvertvirtest,scdim,vcnewdim,vcnewdimmax,vcolddim
692  INTEGER, DIMENSION(:), ALLOCATABLE :: ivgrecvsorted,ivgsendsorted, &
693  ivgsharedsorted
694  TYPE(t_border), POINTER :: pborder,pborderrelated
695  TYPE(t_grid), POINTER :: pgrid
696  TYPE(t_global), POINTER :: global
697  TYPE(t_patch), POINTER :: ppatchrelated,ppatch2,ppatch3
698 
699 ! ******************************************************************************
700 ! Start
701 ! ******************************************************************************
702 
703  global => pregion%global
704 
705  CALL registerfunction(global,'RFLU_SYPE_AugmentCellLists',&
706  'RFLU_ModSymmetryPeriodic.F90')
707 
708  IF ( global%myProcid == masterproc .AND. &
709  global%verbLevel > verbose_none ) THEN
710  WRITE(stdout,'(A,1X,A)') solver_name, &
711  'Adding s/p virtual cells to data structure...'
712  END IF ! global%myProcid
713 
714 ! ******************************************************************************
715 ! Set pointers and initialize variables
716 ! ******************************************************************************
717 
718  pgrid => pregion%grid
719 
720  nvertvirtest = 8*ncellsvirt
721 
722  pborder => pgrid%borders(ppatch%iBorder)
723 
724  ipatchrelated = ppatch%iPatchRelated
725  ppatchrelated => pregion%patches(ipatchrelated)
726 
727  IF ( ppatch%bcType == bc_symmetry ) THEN
728  pborderrelated => pgrid%borders(ppatch%iBorder)
729  ELSE IF ( ppatch%bcType == bc_periodic ) THEN
730  pborderrelated => pgrid%borders(ppatchrelated%iBorder)
731  ELSE
732  CALL errorstop(global,err_reached_default,__line__)
733  END IF ! pPatch%bcType
734 
735  pborder%iBorder = ppatchrelated%iBorder
736  pborder%iRegionGlobal = pregion%iRegionGlobal
737 
738  pborderrelated%nCellsSend = ncellsvirt
739  pborderrelated%nVertSend = 0
740 
741  pborder%nCellsRecv = ncellsvirt
742  pborder%nVertRecv = 0
743  pborder%nVertShared = 0
744 
745  ALLOCATE(pborderrelated%icgSend(pborderrelated%nCellsSend),stat=errorflag)
746  global%error = errorflag
747  IF ( global%error /= err_none ) THEN
748  CALL errorstop(global,err_allocate,__line__,'pBorderRelated%icgSend')
749  END IF ! global%error
750 
751  ALLOCATE(pborderrelated%ivgSend(nvertvirtest),stat=errorflag)
752  global%error = errorflag
753  IF ( global%error /= err_none ) THEN
754  CALL errorstop(global,err_allocate,__line__,'pBorderRelated%ivgSend')
755  END IF ! global%error
756 
757  ALLOCATE(pborder%icgRecv(pborder%nCellsRecv),stat=errorflag)
758  global%error = errorflag
759  IF ( global%error /= err_none ) THEN
760  CALL errorstop(global,err_allocate,__line__,'pBorder%icgRecv')
761  END IF ! global%error
762 
763  ALLOCATE(pborder%ivgRecv(nvertvirtest),stat=errorflag)
764  global%error = errorflag
765  IF ( global%error /= err_none ) THEN
766  CALL errorstop(global,err_allocate,__line__,'pBorder%ivgRecv')
767  END IF ! global%error
768 
769  ALLOCATE(pborder%ivgShared(nvertvirtest),stat=errorflag)
770  global%error = errorflag
771  IF ( global%error /= err_none ) THEN
772  CALL errorstop(global,err_allocate,__line__,'pBorder%ivgShared')
773  END IF ! global%error
774 
775  ALLOCATE(ivgsendsorted(nvertvirtest),stat=errorflag)
776  global%error = errorflag
777  IF ( global%error /= err_none ) THEN
778  CALL errorstop(global,err_allocate,__line__,'ivgSendSorted')
779  END IF ! global%error
780 
781  ALLOCATE(ivgrecvsorted(nvertvirtest),stat=errorflag)
782  global%error = errorflag
783  IF ( global%error /= err_none ) THEN
784  CALL errorstop(global,err_allocate,__line__,'ivgRecvSorted')
785  END IF ! global%error
786 
787  ALLOCATE(ivgsharedsorted(nvertvirtest),stat=errorflag)
788  global%error = errorflag
789  IF ( global%error /= err_none ) THEN
790  CALL errorstop(global,err_allocate,__line__,'ivgSharedSorted')
791  END IF ! global%error
792 
793 ! ******************************************************************************
794 ! Loop over list of virtual cells, add to connectivity lists according to cell
795 ! type
796 ! ******************************************************************************
797 
798  DO icv = 1,ncellsvirt
799  icg = vc(icv)
800 
801  ict = pgrid%cellGlob2Loc(1,icg)
802  icl = pgrid%cellGlob2Loc(2,icg)
803 
804 ! ==============================================================================
805 ! Specify connectivity and set local-to-global mapping
806 ! ==============================================================================
807 
808  SELECT CASE ( ict )
809 
810 ! ------------------------------------------------------------------------------
811 ! Tetrahedra
812 ! ------------------------------------------------------------------------------
813 
814  CASE ( cell_type_tet )
815  pgrid%nCellsTot = pgrid%nCellsTot + 1
816  pgrid%nTetsTot = pgrid%nTetsTot + 1
817 
818  IF ( pgrid%nTetsTot <= pgrid%nTetsMax ) THEN
819  CALL rflu_sype_augmentvertexlists(pregion,pgrid,ppatch,pborder, &
820  pborderrelated,pgrid%tet2v,icl, &
821  pgrid%nTetsTot,ivgsendsorted, &
822  ivgrecvsorted,ivgsharedsorted)
823 
824  pgrid%cellGlob2Loc(1,pgrid%nCellsTot) = cell_type_tet
825  pgrid%cellGlob2Loc(2,pgrid%nCellsTot) = pgrid%nTetsTot
826 
827  pgrid%tet2CellGlob(pgrid%nTetsTot) = pgrid%nCellsTot
828 
829  IF ( ppatch%bcType == bc_symmetry ) THEN
830  CALL swapintegers(pgrid%tet2v(1,pgrid%nTetsTot), &
831  pgrid%tet2v(2,pgrid%nTetsTot))
832  END IF ! pPatch%bcType
833  ELSE
834  CALL errorstop(global,err_exceed_dimens,__line__,'pGrid%tet2v')
835  END IF ! pGrid%nTetsTot
836 
837 ! ------------------------------------------------------------------------------
838 ! Hexahedra
839 ! ------------------------------------------------------------------------------
840 
841  CASE ( cell_type_hex )
842  pgrid%nCellsTot = pgrid%nCellsTot + 1
843  pgrid%nHexsTot = pgrid%nHexsTot + 1
844 
845  IF ( pgrid%nHexsTot <= pgrid%nHexsMax ) THEN
846  CALL rflu_sype_augmentvertexlists(pregion,pgrid,ppatch,pborder, &
847  pborderrelated,pgrid%hex2v,icl, &
848  pgrid%nHexsTot,ivgsendsorted, &
849  ivgrecvsorted,ivgsharedsorted)
850 
851  pgrid%cellGlob2Loc(1,pgrid%nCellsTot) = cell_type_hex
852  pgrid%cellGlob2Loc(2,pgrid%nCellsTot) = pgrid%nHexsTot
853 
854  pgrid%hex2CellGlob(pgrid%nHexsTot) = pgrid%nCellsTot
855 
856  IF ( ppatch%bcType == bc_symmetry ) THEN
857  CALL swapintegers(pgrid%hex2v(2,pgrid%nHexsTot), &
858  pgrid%hex2v(4,pgrid%nHexsTot))
859  CALL swapintegers(pgrid%hex2v(6,pgrid%nHexsTot), &
860  pgrid%hex2v(8,pgrid%nHexsTot))
861  END IF ! pPatch%bcType
862  ELSE
863  CALL errorstop(global,err_exceed_dimens,__line__,'pGrid%hex2v')
864  END IF ! pGrid%nHexsTot
865 
866 ! ------------------------------------------------------------------------------
867 ! Prisms
868 ! ------------------------------------------------------------------------------
869 
870  CASE ( cell_type_pri )
871  pgrid%nCellsTot = pgrid%nCellsTot + 1
872  pgrid%nPrisTot = pgrid%nPrisTot + 1
873 
874  IF ( pgrid%nPrisTot <= pgrid%nPrisMax ) THEN
875  CALL rflu_sype_augmentvertexlists(pregion,pgrid,ppatch,pborder, &
876  pborderrelated,pgrid%pri2v,icl, &
877  pgrid%nPrisTot,ivgsendsorted, &
878  ivgrecvsorted,ivgsharedsorted)
879 
880  pgrid%cellGlob2Loc(1,pgrid%nCellsTot) = cell_type_pri
881  pgrid%cellGlob2Loc(2,pgrid%nCellsTot) = pgrid%nPrisTot
882 
883  pgrid%pri2CellGlob(pgrid%nPrisTot) = pgrid%nCellsTot
884 
885  IF ( ppatch%bcType == bc_symmetry ) THEN
886  CALL swapintegers(pgrid%pri2v(2,pgrid%nPrisTot), &
887  pgrid%pri2v(3,pgrid%nPrisTot))
888  CALL swapintegers(pgrid%pri2v(5,pgrid%nPrisTot), &
889  pgrid%pri2v(6,pgrid%nPrisTot))
890  END IF ! pPatch%bcType
891  ELSE
892  CALL errorstop(global,err_exceed_dimens,__line__,'pGrid%pri2v')
893  END IF ! pGrid%nPrisTot
894 
895 ! ------------------------------------------------------------------------------
896 ! Pyramids
897 ! ------------------------------------------------------------------------------
898 
899  CASE ( cell_type_pyr )
900  pgrid%nCellsTot = pgrid%nCellsTot + 1
901  pgrid%nPyrsTot = pgrid%nPyrsTot + 1
902 
903  IF ( pgrid%nPyrsTot <= pgrid%nPyrsMax ) THEN
904  CALL rflu_sype_augmentvertexlists(pregion,pgrid,ppatch,pborder, &
905  pborderrelated,pgrid%pyr2v,icl, &
906  pgrid%nPyrsTot,ivgsendsorted, &
907  ivgrecvsorted,ivgsharedsorted)
908 
909  pgrid%cellGlob2Loc(1,pgrid%nCellsTot) = cell_type_pyr
910  pgrid%cellGlob2Loc(2,pgrid%nCellsTot) = pgrid%nPyrsTot
911 
912  pgrid%pyr2CellGlob(pgrid%nPyrsTot) = pgrid%nCellsTot
913 
914  IF ( ppatch%bcType == bc_symmetry ) THEN
915  CALL swapintegers(pgrid%pyr2v(2,pgrid%nPyrsTot), &
916  pgrid%pyr2v(4,pgrid%nPyrsTot))
917  END IF ! pPatch%bcType
918  ELSE
919  CALL errorstop(global,err_exceed_dimens,__line__,'pGrid%pyr2v')
920  END IF ! pGrid%nPyrsTot
921 
922 ! ------------------------------------------------------------------------------
923 ! Default
924 ! ------------------------------------------------------------------------------
925 
926  CASE default
927  CALL errorstop(global,err_reached_default,__line__)
928  END SELECT ! pGridSerial%cellGlob2Loc
929 
930 ! ==============================================================================
931 ! Build cell communication lists
932 ! ==============================================================================
933 
934  pborderrelated%icgSend(icv) = icg
935 
936  pborder%icgRecv(icv) = pgrid%nCellsTot
937  END DO ! icl
938 
939 ! ******************************************************************************
940 ! Add faces of virtual cells to patch face lists
941 ! ******************************************************************************
942 
943  DO ipatch2 = 1,pgrid%nPatches
944  ppatch2 => pregion%patches(ipatch2)
945 
946  ALLOCATE(ppatch2%bf2cSorted(ppatch2%nBFaces),stat=errorflag)
947  global%error = errorflag
948  IF ( global%error /= err_none ) THEN
949  CALL errorstop(global,err_allocate,__line__,'pPatch%bf2cSorted')
950  END IF ! global%error
951 
952  ALLOCATE(ppatch2%bf2cSortedKeys(ppatch2%nBFaces),stat=errorflag)
953  global%error = errorflag
954  IF ( global%error /= err_none ) THEN
955  CALL errorstop(global,err_allocate,__line__,'pPatch%bf2cSortedKeys')
956  END IF ! global%error
957 
958  ALLOCATE(ppatch2%bf2vSorted(4,ppatch2%nBFaces),stat=errorflag)
959  global%error = errorflag
960  IF ( global%error /= err_none ) THEN
961  CALL errorstop(global,err_allocate,__line__,'pPatch%bf2vSorted')
962  END IF ! global%error
963 
964  DO ifl = 1,ppatch2%nBFaces
965  ppatch2%bf2cSorted(ifl) = ppatch2%bf2c(ifl)
966  ppatch2%bf2cSortedKeys(ifl) = ifl
967  END DO ! ifl
968 
969  CALL quicksortintegerinteger(ppatch2%bf2cSorted,ppatch2%bf2cSortedKeys, &
970  ppatch2%nBFaces)
971 
972  DO ifl = 1,ppatch2%nBFaces
973  ifl2 = ppatch2%bf2cSortedKeys(ifl)
974 
975  ppatch2%bf2vSorted(1,ifl) = ppatch2%bf2v(1,ifl2)
976  ppatch2%bf2vSorted(2,ifl) = ppatch2%bf2v(2,ifl2)
977  ppatch2%bf2vSorted(3,ifl) = ppatch2%bf2v(3,ifl2)
978  ppatch2%bf2vSorted(4,ifl) = ppatch2%bf2v(4,ifl2)
979  END DO ! ifl
980  END DO ! iPatch2
981 
982 ! ******************************************************************************
983 ! Loop over cells to be sent
984 ! ******************************************************************************
985 
986  DO icl = 1,pborderrelated%nCellsSend
987  icg = pborderrelated%icgSend(icl)
988 
989 ! ==============================================================================
990 ! Loop over patches
991 ! ==============================================================================
992 
993  DO ipatch2 = 1,pgrid%nPatches
994  ppatch2 => pregion%patches(ipatch2)
995 
996 ! ------------------------------------------------------------------------------
997 ! If current patch is not the same as patch associated with symmetry or
998 ! periodic patch, find out whether cell to be sent is adjacent to current
999 ! patch
1000 ! ------------------------------------------------------------------------------
1001 
1002  IF ( (ppatch2%iPatchGlobal /= ppatch%iPatchGlobal) .AND. &
1003  (ppatch2%nbMap(ppatch%iPatchGlobal) .EQV. .true. ) ) THEN
1004  CALL binarysearchinteger(ppatch2%bf2cSorted(1:ppatch2%nBFaces), &
1005  ppatch2%nBFaces,icg,iloc)
1006 
1007 ! ----- Cell found, so adjacent to current patch -------------------------------
1008 
1009  IF ( iloc /= element_not_found ) THEN
1010 
1011 ! ------- Triangle face
1012 
1013  IF ( ppatch2%bf2vSorted(4,iloc) == vert_none ) THEN
1014  IF ( ppatch2%nBTrisTot < ppatch2%nBTrisMax ) THEN
1015  ppatch2%nBTrisTot = ppatch2%nBTrisTot + 1
1016  ppatch2%nBFacesTot = ppatch2%nBFacesTot + 1
1017  ELSE
1018  CALL errorstop(global,err_exceed_dimens,__line__,'tri2v')
1019  END IF ! pPatch2%nBTrisTot
1020 
1021  DO ivl2 = 1,3
1022  ivl = ppatch2%bf2vSorted(ivl2,iloc)
1023  ivg = ppatch2%bv(ivl)
1024 
1025  CALL binarysearchinteger(ivgsendsorted(1:pborderrelated%nVertSend), &
1026  pborderrelated%nVertSend,ivg,iloc2)
1027 
1028  IF ( iloc2 /= element_not_found ) THEN
1029  ppatch2%bTri2v(ivl2,ppatch2%nBTrisTot) = ivgrecvsorted(iloc2)
1030  ELSE
1031  IF ( ppatch%bcType == bc_symmetry ) THEN
1032  CALL binarysearchinteger(ivgsharedsorted(1:pborder%nVertShared), &
1033  pborder%nVertShared,ivg,iloc2)
1034  ELSE IF ( ppatch%bcType == bc_periodic ) THEN
1035  ipatchrelated = ppatch%iPatchRelated
1036 
1037  ppatch3 => pregion%patches(ipatchrelated)
1038 
1039  CALL binarysearchinteger(ppatch3%bv(1:ppatch3%nBVert), &
1040  ppatch3%nBVert,ivg,iloc2)
1041  ELSE
1042  CALL errorstop(global,err_reached_default,__line__)
1043  END IF ! pPatch%bcType
1044 
1045  IF ( iloc2 /= element_not_found ) THEN
1046  IF ( ppatch%bcType == bc_symmetry ) THEN
1047  ppatch2%bTri2v(ivl2,ppatch2%nBTrisTot) = ivg
1048  ELSE IF ( ppatch%bcType == bc_periodic ) THEN
1049  ivlm = ppatch3%bv2bv(iloc2)
1050  ivgm = ppatch%bv(ivlm)
1051 
1052  ppatch2%bTri2v(ivl2,ppatch2%nBTrisTot) = ivgm
1053  END IF ! pPatch%bcType
1054  ELSE
1055  CALL errorstop(global,err_vertex_not_found,__line__)
1056  END IF ! iLoc2
1057  END IF ! iLoc2
1058  END DO ! ivl2
1059 
1060  IF ( ppatch%bcType == bc_symmetry ) THEN
1061  CALL swapintegers(ppatch2%bTri2v(2,ppatch2%nBTrisTot), &
1062  ppatch2%bTri2v(3,ppatch2%nBTrisTot))
1063  END IF ! pPatch%bcType
1064 
1065 ! ------- Quadrilateral face
1066 
1067  ELSE
1068  IF ( ppatch2%nBQuadsTot < ppatch2%nBQuadsMax ) THEN
1069  ppatch2%nBQuadsTot = ppatch2%nBQuadsTot + 1
1070  ppatch2%nBFacesTot = ppatch2%nBFacesTot + 1
1071  ELSE
1072  CALL errorstop(global,err_exceed_dimens,__line__,'quad2v')
1073  END IF ! pPatch2%nBQuadsTot
1074 
1075  DO ivl2 = 1,4
1076  ivl = ppatch2%bf2vSorted(ivl2,iloc)
1077  ivg = ppatch2%bv(ivl)
1078 
1079  CALL binarysearchinteger(ivgsendsorted(1:pborderrelated%nVertSend), &
1080  pborderrelated%nVertSend,ivg,iloc2)
1081 
1082  IF ( iloc2 /= element_not_found ) THEN
1083  ppatch2%bQuad2v(ivl2,ppatch2%nBQuadsTot) = ivgrecvsorted(iloc2)
1084  ELSE
1085  IF ( ppatch%bcType == bc_symmetry ) THEN
1086  CALL binarysearchinteger(ivgsharedsorted(1:pborder%nVertShared), &
1087  pborder%nVertShared,ivg,iloc2)
1088  ELSE IF ( ppatch%bcType == bc_periodic ) THEN
1089  ipatchrelated = ppatch%iPatchRelated
1090 
1091  ppatch3 => pregion%patches(ipatchrelated)
1092 
1093  CALL binarysearchinteger(ppatch3%bv(1:ppatch3%nBVert), &
1094  ppatch3%nBVert,ivg,iloc2)
1095  ELSE
1096  CALL errorstop(global,err_reached_default,__line__)
1097  END IF ! pPatch%bcType
1098 
1099  IF ( iloc2 /= element_not_found ) THEN
1100  IF ( ppatch%bcType == bc_symmetry ) THEN
1101  ppatch2%bQuad2v(ivl2,ppatch2%nBQuadsTot) = ivg
1102  ELSE IF ( ppatch%bcType == bc_periodic ) THEN
1103  ivlm = ppatch3%bv2bv(iloc2)
1104  ivgm = ppatch%bv(ivlm)
1105 
1106  ppatch2%bQuad2v(ivl2,ppatch2%nBQuadsTot) = ivgm
1107  END IF ! pPatch%bcType
1108  ELSE
1109  CALL errorstop(global,err_vertex_not_found,__line__)
1110  END IF ! iLoc2
1111  END IF ! iLoc2
1112  END DO ! ivl2
1113 
1114  IF ( ppatch%bcType == bc_symmetry ) THEN
1115  CALL swapintegers(ppatch2%bQuad2v(2,ppatch2%nBQuadsTot), &
1116  ppatch2%bQuad2v(4,ppatch2%nBQuadsTot))
1117  END IF ! pPatch%bcType
1118  END IF ! pPatch2%bf2v
1119  END IF ! iLoc
1120  END IF ! pPatch2%iPatchGlobal
1121  END DO ! iPatch
1122  END DO ! icl
1123 
1124  DO ipatch2 = 1,pgrid%nPatches
1125  ppatch2 => pregion%patches(ipatch2)
1126 
1127  DEALLOCATE(ppatch2%bf2cSorted,stat=errorflag)
1128  global%error = errorflag
1129  IF ( global%error /= err_none ) THEN
1130  CALL errorstop(global,err_deallocate,__line__,'pPatch%bf2cSorted')
1131  END IF ! global%error
1132  END DO ! iPatch2
1133 
1134 ! ******************************************************************************
1135 ! Deallocate temporary memory
1136 ! ******************************************************************************
1137 
1138  DEALLOCATE(ivgsendsorted,stat=errorflag)
1139  global%error = errorflag
1140  IF ( global%error /= err_none ) THEN
1141  CALL errorstop(global,err_deallocate,__line__,'ivgSendSorted')
1142  END IF ! global%error
1143 
1144  DEALLOCATE(ivgrecvsorted,stat=errorflag)
1145  global%error = errorflag
1146  IF ( global%error /= err_none ) THEN
1147  CALL errorstop(global,err_deallocate,__line__,'ivgRecvSorted')
1148  END IF ! global%error
1149 
1150 ! ******************************************************************************
1151 ! End
1152 ! ******************************************************************************
1153 
1154  IF ( global%myProcid == masterproc .AND. &
1155  global%verbLevel > verbose_none ) THEN
1156  WRITE(stdout,'(A,1X,A)') solver_name, &
1157  'Adding s/p virtual cells to data structure done.'
1158  END IF ! global%myProcid
1159 
1160  CALL deregisterfunction(global)
1161 
1162 END SUBROUTINE rflu_sype_augmentcelllists
1163 
1164 
1165 
1166 
1167 
1168 
1169 ! ******************************************************************************
1170 !
1171 ! Purpose: Loop over vertices of given cell and update lists relating to
1172 ! vertices.
1173 !
1174 ! Description: None.
1175 !
1176 ! Input:
1177 ! pRegion Pointer to region
1178 ! pGrid Pointer to grid
1179 ! pPatch Pointer to patch
1180 ! pBorder Pointer to border
1181 ! pBorderRelated Pointer to related border
1182 ! x2v Cell-to-vertex connectivity
1183 ! icl Local cell index
1184 ! nXTot Total number of cells of given type
1185 ! ivgSendSorted List of vertices to be sent for given border
1186 ! ivgRecvSorted List of vertices to be received for given border
1187 ! ivgSharedSorted List of vertices to be shared for given border
1188 !
1189 ! Output:
1190 ! x2v Cell-to-vertex connectivity
1191 ! nXTot Total number of cells of given type
1192 ! ivgSendSorted List of vertices to be sent for given border
1193 ! ivgRecvSorted List of vertices to be received for given border
1194 ! ivgSharedSorted List of vertices to be shared for given border
1195 !
1196 ! Notes: None.
1197 !
1198 ! ******************************************************************************
1199 
1200 SUBROUTINE rflu_sype_augmentvertexlists(pRegion,pGrid,pPatch,pBorder, &
1201  pborderrelated,x2v,icl,nxtot, &
1202  ivgsendsorted,ivgrecvsorted, &
1203  ivgsharedsorted)
1204 
1205  IMPLICIT NONE
1206 
1207 ! ******************************************************************************
1208 ! Declarations and definitions
1209 ! ******************************************************************************
1210 
1211 ! ==============================================================================
1212 ! Arguments
1213 ! ==============================================================================
1214 
1215  INTEGER, INTENT(IN) :: icl
1216  INTEGER, INTENT(INOUT) :: nxtot
1217  INTEGER, DIMENSION(:) :: ivgsendsorted,ivgsharedsorted,ivgrecvsorted
1218  INTEGER, DIMENSION(:,:), POINTER :: x2v
1219  REAL(RFREAL), DIMENSION(XCOORD:XYZMAG) :: xyz,xyzt
1220  TYPE(t_border), POINTER :: pborder,pborderrelated
1221  TYPE(t_grid), POINTER :: pgrid
1222  TYPE(t_patch), POINTER :: ppatch
1223  TYPE(t_region), POINTER :: pregion
1224 
1225 ! ==============================================================================
1226 ! Locals
1227 ! ==============================================================================
1228 
1229  INTEGER :: iloc,iloc2,ipatchrelated,ivg,ivgm,ivl
1230  REAL(RFREAL) :: nxpc,nypc,nzpc,xpc,xvg,ypc,yvg,zpc,zvg
1231  TYPE(t_global), POINTER :: global
1232  TYPE(t_patch), POINTER :: ppatchrelated
1233 
1234 ! ******************************************************************************
1235 ! Start
1236 ! ******************************************************************************
1237 
1238  global => pregion%global
1239 
1240  CALL registerfunction(global,'RFLU_SYPE_AugmentVertexLists',&
1241  'RFLU_ModSymmetryPeriodic.F90')
1242 
1243 ! ******************************************************************************
1244 ! Get patch geometry and set patch index and pointer
1245 ! ******************************************************************************
1246 
1247  xpc = ppatch%pc(xcoord)
1248  ypc = ppatch%pc(ycoord)
1249  zpc = ppatch%pc(zcoord)
1250 
1251  nxpc = ppatch%pn(xcoord)
1252  nypc = ppatch%pn(ycoord)
1253  nzpc = ppatch%pn(zcoord)
1254 
1255  IF ( ppatch%bcType == bc_periodic ) THEN
1256  ipatchrelated = ppatch%iPatchRelated
1257 
1258  ppatchrelated => pregion%patches(ipatchrelated)
1259  ELSE
1260  ipatchrelated = crazy_value_int
1261 
1262  nullify(ppatchrelated)
1263  END IF ! pPatch%bcType
1264 
1265 ! ******************************************************************************
1266 ! Loop over vertices of cell
1267 ! ******************************************************************************
1268 
1269  DO ivl = 1,SIZE(x2v,1)
1270  ivg = x2v(ivl,icl)
1271 
1272 ! ==============================================================================
1273 ! Check whether vertex is on patch for symmetry patches or on related patch
1274 ! for periodic patches
1275 ! ==============================================================================
1276 
1277  IF ( ppatch%bcType == bc_symmetry ) THEN
1278  CALL binarysearchinteger(ppatch%bv(1:ppatch%nBVert),ppatch%nBVert,ivg, &
1279  iloc)
1280  ELSE IF ( ppatch%bcType == bc_periodic ) THEN
1281  CALL binarysearchinteger(ppatchrelated%bv(1:ppatchrelated%nBVert), &
1282  ppatchrelated%nBVert,ivg,iloc)
1283  ELSE
1284  CALL errorstop(global,err_reached_default,__line__)
1285  END IF ! pPatch%bcType
1286 
1287 ! ==============================================================================
1288 ! Vertex not on patch
1289 ! ==============================================================================
1290 
1291  IF ( iloc == element_not_found ) THEN
1292 
1293 ! ------------------------------------------------------------------------------
1294 ! Check whether already in list of vertices to be sent
1295 ! ------------------------------------------------------------------------------
1296 
1297  IF ( pborderrelated%nVertSend > 1 ) THEN
1298  CALL binarysearchinteger(ivgsendsorted(1:pborderrelated%nVertSend), &
1299  pborderrelated%nVertSend,ivg,iloc2)
1300  ELSE
1301  iloc2 = element_not_found
1302  END IF ! pBorderRelated%nVertSend
1303 
1304 ! ------------------------------------------------------------------------------
1305 ! Vertex not in list of vertices to be sent, so add. Compute new coordinates
1306 ! and add vertex to lists
1307 ! ------------------------------------------------------------------------------
1308 
1309  IF ( iloc2 == element_not_found ) THEN
1310  IF ( pgrid%nVertTot < pgrid%nVertMax ) THEN
1311  pgrid%nVertTot = pgrid%nVertTot + 1
1312  x2v(ivl,nxtot) = pgrid%nVertTot
1313  ELSE
1314  CALL errorstop(global,err_exceed_dimens,__line__,'x2v')
1315  END IF ! pGrid%nVertTot
1316 
1317  IF ( ppatch%bcType == bc_symmetry ) THEN
1318  xvg = pgrid%xyz(xcoord,ivg)
1319  yvg = pgrid%xyz(ycoord,ivg)
1320  zvg = pgrid%xyz(zcoord,ivg)
1321 
1322  CALL reflectposition(nxpc,nypc,nzpc,xpc,ypc,zpc,xvg,yvg,zvg)
1323 
1324  pgrid%xyz(xcoord,pgrid%nVertTot) = xvg
1325  pgrid%xyz(ycoord,pgrid%nVertTot) = yvg
1326  pgrid%xyz(zcoord,pgrid%nVertTot) = zvg
1327  ELSE IF ( ppatch%bcType == bc_periodic ) THEN
1328  xyz(xcoord) = pgrid%xyz(xcoord,ivg)
1329  xyz(ycoord) = pgrid%xyz(ycoord,ivg)
1330  xyz(zcoord) = pgrid%xyz(zcoord,ivg)
1331  xyz(xyzmag) = 1.0_rfreal
1332 
1333  xyzt = matmul(ppatchrelated%tm,xyz)
1334 
1335  pgrid%xyz(xcoord,pgrid%nVertTot) = xyzt(xcoord)
1336  pgrid%xyz(ycoord,pgrid%nVertTot) = xyzt(ycoord)
1337  pgrid%xyz(zcoord,pgrid%nVertTot) = xyzt(zcoord)
1338  ELSE
1339  CALL errorstop(global,err_reached_default,__line__)
1340  END IF ! pPatch%bcType
1341 
1342  pborderrelated%nVertSend = pborderrelated%nVertSend + 1
1343 
1344  pborderrelated%ivgSend(pborderrelated%nVertSend) = ivg
1345  ivgsendsorted(pborderrelated%nVertSend) = ivg
1346 
1347  pborder%nVertRecv = pborder%nVertRecv + 1
1348 
1349  pborder%ivgRecv(pborder%nVertRecv) = pgrid%nVertTot
1350  ivgrecvsorted(pborder%nVertRecv) = pgrid%nVertTot
1351 
1352  IF ( pborder%nVertRecv > 1 ) THEN
1353  CALL quicksortintegerinteger(ivgsendsorted(1:pborder%nVertRecv), &
1354  ivgrecvsorted(1:pborder%nVertRecv), &
1355  pborder%nVertRecv)
1356  END IF ! pBorder%nVertRecv
1357 
1358 ! ------------------------------------------------------------------------------
1359 ! Vertex already in list of vertices to be sent
1360 ! ------------------------------------------------------------------------------
1361 
1362  ELSE
1363  x2v(ivl,nxtot) = ivgrecvsorted(iloc2)
1364  END IF ! iLoc2
1365 
1366 ! ==============================================================================
1367 ! Vertex on patch
1368 ! ==============================================================================
1369 
1370  ELSE
1371  IF ( ppatch%bcType == bc_symmetry ) THEN
1372  x2v(ivl,nxtot) = ivg
1373  ELSE IF ( ppatch%bcType == bc_periodic ) THEN
1374  CALL rflu_sype_imposevertexmap(pregion,ppatchrelated,ivg,ivgm)
1375 
1376  x2v(ivl,nxtot) = ivgm
1377  ELSE
1378  CALL errorstop(global,err_reached_default,__line__)
1379  END IF ! pPatch%bcType
1380 
1381 ! ------------------------------------------------------------------------------
1382 ! Check whether already in list of vertices to be shared
1383 ! ------------------------------------------------------------------------------
1384 
1385  IF ( pborder%nVertShared > 1 ) THEN
1386  IF ( ppatch%bcType == bc_symmetry ) THEN
1387  CALL binarysearchinteger(ivgsharedsorted(1:pborder%nVertShared), &
1388  pborder%nVertShared,ivg,iloc2)
1389  ELSE IF ( ppatch%bcType == bc_periodic ) THEN
1390  CALL binarysearchinteger(ivgsharedsorted(1:pborder%nVertShared), &
1391  pborder%nVertShared,ivgm,iloc2)
1392  ELSE
1393  CALL errorstop(global,err_reached_default,__line__)
1394  END IF ! pPatch%bcType
1395  ELSE
1396  iloc2 = element_not_found
1397  END IF ! pBorder%nVertShared
1398 
1399 ! ------------------------------------------------------------------------------
1400 ! Vertex not in list of vertices to be shared, so add
1401 ! ------------------------------------------------------------------------------
1402 
1403  IF ( iloc2 == element_not_found ) THEN
1404  pborder%nVertShared = pborder%nVertShared + 1
1405 
1406  IF ( ppatch%bcType == bc_symmetry ) THEN
1407  pborder%ivgShared(pborder%nVertShared) = ivg
1408  ivgsharedsorted(pborder%nVertShared) = ivg
1409  ELSE IF ( ppatch%bcType == bc_periodic ) THEN
1410  pborder%ivgShared(pborder%nVertShared) = ivgm
1411  ivgsharedsorted(pborder%nVertShared) = ivgm
1412  ELSE
1413  CALL errorstop(global,err_reached_default,__line__)
1414  END IF ! pPatch%bcType
1415 
1416  IF ( pborder%nVertShared > 1 ) THEN
1417  CALL quicksortinteger(ivgsharedsorted(1:pborder%nVertShared), &
1418  pborder%nVertShared)
1419  END IF ! pBorder%nVertShared
1420  END IF ! iLoc2
1421  END IF ! iLoc
1422  END DO ! ivl
1423 
1424 ! ******************************************************************************
1425 ! End
1426 ! ******************************************************************************
1427 
1428  CALL deregisterfunction(global)
1429 
1430 END SUBROUTINE rflu_sype_augmentvertexlists
1431 
1432 
1433 
1434 
1435 
1436 
1437 
1438 ! ******************************************************************************
1439 !
1440 ! Purpose: Build list of virtual cells adjacent to periodic and symmetry patches
1441 ! on a parallel region.
1442 !
1443 ! Description: None.
1444 !
1445 ! Input:
1446 ! pRegion Pointer to parallel region
1447 ! pRegionSerial Pointer to serial region
1448 !
1449 ! Output: None.
1450 !
1451 ! Notes: None.
1452 !
1453 ! ******************************************************************************
1454 
1455  SUBROUTINE rflu_sype_buildp2vclist(pRegion,pRegionSerial)
1456 
1457  IMPLICIT NONE
1458 
1459 ! ******************************************************************************
1460 ! Declarations and definitions
1461 ! ******************************************************************************
1462 
1463 ! ==============================================================================
1464 ! Arguments
1465 ! ==============================================================================
1466 
1467  TYPE(t_region), POINTER :: pregion,pregionserial
1468 
1469 ! ==============================================================================
1470 ! Locals
1471 ! ==============================================================================
1472 
1473  INTEGER :: errorflag,iborder,icg,icgs,icl,iloc,ipatch,ipatchserial
1474  INTEGER, DIMENSION(:), ALLOCATABLE :: ipatchglob2ipatchloc
1475  TYPE(t_global), POINTER :: global
1476  TYPE(t_border), POINTER :: pborder
1477  TYPE(t_grid), POINTER :: pgrid,pgridserial
1478  TYPE(t_patch), POINTER :: ppatch,ppatchserial
1479 
1480 ! ******************************************************************************
1481 ! Start
1482 ! ******************************************************************************
1483 
1484  global => pregion%global
1485 
1486  CALL registerfunction(global,'RFLU_SYPE_BuildP2VCList',&
1487  'RFLU_ModSymmetryPeriodic.F90')
1488 
1489  IF ( global%myProcid == masterproc .AND. &
1490  global%verbLevel > verbose_low ) THEN
1491  WRITE(stdout,'(A,1X,A)') solver_name, &
1492  'Building patch-to-virtual-cell list...'
1493  WRITE(stdout,'(A,3X,A,1X,I5.5)') solver_name,'Global region:', &
1494  pregion%iRegionGlobal
1495  END IF ! global%verbLevel
1496 
1497 ! ******************************************************************************
1498 ! Check and set grid pointer
1499 ! ******************************************************************************
1500 
1501  IF ( pregion%iRegionGlobal == 0 ) THEN ! Defensive coding
1502  CALL errorstop(global,err_region_id_invalid,__line__)
1503  END IF ! pRegion%iRegionGlobal
1504 
1505  pgrid => pregion%grid
1506  pgridserial => pregionserial%grid
1507 
1508 ! ******************************************************************************
1509 ! Allocate temporary memory for and build mapping of global to local patch
1510 ! indices
1511 ! ******************************************************************************
1512 
1513  ALLOCATE(ipatchglob2ipatchloc(pgridserial%nPatches),stat=errorflag)
1514  global%error = errorflag
1515  IF ( global%error /= err_none ) THEN
1516  CALL errorstop(global,err_allocate,__line__,'iPatchGlob2iPatchLoc')
1517  END IF ! global%error
1518 
1519  DO ipatchserial = 1,pgridserial%nPatches
1520  ipatchglob2ipatchloc(ipatchserial) = crazy_value_int
1521  END DO ! iPatchSerial
1522 
1523  DO ipatch = 1,pgrid%nPatches
1524  ppatch => pregion%patches(ipatch)
1525 
1526  ipatchglob2ipatchloc(ppatch%iPatchGlobal) = ipatch
1527  END DO ! iPatch
1528 
1529 ! ******************************************************************************
1530 ! Initialize number of virtual cells for each patch. NOTE should not be
1531 ! necessary, but done for safety
1532 ! ******************************************************************************
1533 
1534  DO ipatch = 1,pgrid%nPatches
1535  ppatch => pregion%patches(ipatch)
1536 
1537  ppatch%nBCellsVirt = 0
1538  END DO ! iPatch
1539 
1540 ! ******************************************************************************
1541 ! Determine number of virtual cells for periodic and symmetry patches
1542 ! ******************************************************************************
1543 
1544  DO iborder = 1,pgrid%nBorders
1545  pborder => pgrid%borders(iborder)
1546 
1547  DO icl = 1,pborder%nCellsRecv
1548  icg = pborder%icgRecv(icl)
1549  icgs = pgrid%pc2sc(icg)
1550 
1551  patchserialloop: DO ipatchserial = 1,pgridserial%nPatches
1552  ppatchserial => pregionserial%patches(ipatchserial)
1553 
1554  IF ( ppatchserial%bcType == bc_periodic .OR. &
1555  ppatchserial%bcType == bc_symmetry ) THEN
1556  CALL binarysearchinteger(ppatchserial%bvc, &
1557  ppatchserial%nBCellsVirt,icgs,iloc)
1558 
1559  IF ( iloc /= element_not_found ) THEN
1560  ipatch = ipatchglob2ipatchloc(ipatchserial)
1561  ppatch => pregion%patches(ipatch)
1562 
1563  ppatch%nBCellsVirt = ppatch%nBCellsVirt + 1
1564 
1565  EXIT patchserialloop
1566  END IF ! iLoc
1567  END IF ! pPatchSerial%bcType
1568  END DO patchserialloop
1569  END DO ! icl
1570  END DO ! iBorder
1571 
1572 ! ******************************************************************************
1573 ! Allocate memory and set counters back to zero
1574 ! ******************************************************************************
1575 
1576  CALL rflu_sype_createp2vclist(pregion)
1577 
1578  DO ipatch = 1,pgrid%nPatches
1579  ppatch => pregion%patches(ipatch)
1580 
1581  ppatch%nBCellsVirt = 0
1582  END DO ! iPatch
1583 
1584 ! ******************************************************************************
1585 ! Find and store virtual cells for periodic and symmetry patches
1586 ! ******************************************************************************
1587 
1588  DO iborder = 1,pgrid%nBorders
1589  pborder => pgrid%borders(iborder)
1590 
1591  DO icl = 1,pborder%nCellsRecv
1592  icg = pborder%icgRecv(icl)
1593  icgs = pgrid%pc2sc(icg)
1594 
1595  patchserialloop2: DO ipatchserial = 1,pgridserial%nPatches
1596  ppatchserial => pregionserial%patches(ipatchserial)
1597 
1598  IF ( ppatchserial%bcType == bc_periodic .OR. &
1599  ppatchserial%bcType == bc_symmetry ) THEN
1600  CALL binarysearchinteger(ppatchserial%bvc, &
1601  ppatchserial%nBCellsVirt,icgs,iloc)
1602 
1603  IF ( iloc /= element_not_found ) THEN
1604  ipatch = ipatchglob2ipatchloc(ipatchserial)
1605  ppatch => pregion%patches(ipatch)
1606 
1607  ppatch%nBCellsVirt = ppatch%nBCellsVirt + 1
1608  ppatch%bvc(ppatch%nBCellsVirt) = icg
1609 
1610  EXIT patchserialloop2
1611  END IF ! iLoc
1612  END IF ! pPatchSerial%bcType
1613  END DO patchserialloop2
1614  END DO ! icl
1615  END DO ! iBorder
1616 
1617 ! ******************************************************************************
1618 ! Deallocate temporary memory
1619 ! ******************************************************************************
1620 
1621  DEALLOCATE(ipatchglob2ipatchloc,stat=errorflag)
1622  global%error = errorflag
1623  IF ( global%error /= err_none ) THEN
1624  CALL errorstop(global,err_deallocate,__line__,'iPatchGlob2iPatchLoc')
1625  END IF ! global%error
1626 
1627 ! ******************************************************************************
1628 ! End
1629 ! ******************************************************************************
1630 
1631  IF ( global%myProcid == masterproc .AND. &
1632  global%verbLevel > verbose_low ) THEN
1633  WRITE(stdout,'(A,1X,A)') solver_name, &
1634  'Building patch-to-virtual-cell list done.'
1635  END IF ! global%verbLevel
1636 
1637  CALL deregisterfunction(global)
1638 
1639  END SUBROUTINE rflu_sype_buildp2vclist
1640 
1641 
1642 
1643 
1644 
1645 
1646 
1647 ! ******************************************************************************
1648 !
1649 ! Purpose: Build list of virtual cells adjacent to periodic and symmetry patches
1650 ! on serial region.
1651 !
1652 ! Description: None.
1653 !
1654 ! Input:
1655 ! pRegion Pointer to region
1656 !
1657 ! Output: None.
1658 !
1659 ! Notes:
1660 ! 1. This routine can only be used on the serial region.
1661 !
1662 ! ******************************************************************************
1663 
1664  SUBROUTINE rflu_sype_buildp2vclistserial(pRegion)
1665 
1666  IMPLICIT NONE
1667 
1668 ! ******************************************************************************
1669 ! Declarations and definitions
1670 ! ******************************************************************************
1671 
1672 ! ==============================================================================
1673 ! Arguments
1674 ! ==============================================================================
1675 
1676  TYPE(t_region), POINTER :: pregion
1677 
1678 ! ==============================================================================
1679 ! Locals
1680 ! ==============================================================================
1681 
1682  INTEGER :: errorflag,iborder,icl,ipatch
1683  TYPE(t_global), POINTER :: global
1684  TYPE(t_border), POINTER :: pborder
1685  TYPE(t_grid), POINTER :: pgrid
1686  TYPE(t_patch), POINTER :: ppatch
1687 
1688 ! ******************************************************************************
1689 ! Start
1690 ! ******************************************************************************
1691 
1692  global => pregion%global
1693 
1694  CALL registerfunction(global,'RFLU_SYPE_BuildP2VCListSerial',&
1695  'RFLU_ModSymmetryPeriodic.F90')
1696 
1697  IF ( global%myProcid == masterproc .AND. &
1698  global%verbLevel > verbose_low ) THEN
1699  WRITE(stdout,'(A,1X,A)') solver_name, &
1700  'Building patch-to-virtual-cell list...'
1701  WRITE(stdout,'(A,3X,A,1X,I5.5)') solver_name,'Global region:', &
1702  pregion%iRegionGlobal
1703  END IF ! global%verbLevel
1704 
1705 ! ******************************************************************************
1706 ! Check and set grid pointer
1707 ! ******************************************************************************
1708 
1709  IF ( pregion%iRegionGlobal /= 0 ) THEN ! Defensive coding
1710  CALL errorstop(global,err_region_id_invalid,__line__)
1711  END IF ! pRegion%iRegionGlobal
1712 
1713  pgrid => pregion%grid
1714 
1715 ! ******************************************************************************
1716 ! Build list
1717 ! ******************************************************************************
1718 
1719  DO ipatch = 1,pgrid%nPatches
1720  ppatch => pregion%patches(ipatch)
1721 
1722  IF ( ppatch%bcType == bc_periodic .OR. &
1723  ppatch%bcType == bc_symmetry ) THEN
1724  iborder = ppatch%iBorder
1725  pborder => pgrid%borders(iborder)
1726 
1727  ppatch%nBCellsVirt = pborder%nCellsRecv
1728 
1729  ALLOCATE(ppatch%bvc(ppatch%nBCellsVirt),stat=errorflag)
1730  global%error = errorflag
1731  IF ( global%error /= err_none ) THEN
1732  CALL errorstop(global,err_allocate,__line__,'pPatch%bvc')
1733  END IF ! global%error
1734 
1735  DO icl = 1,ppatch%nBCellsVirt
1736  ppatch%bvc(icl) = pborder%icgRecv(icl)
1737  END DO ! icl
1738  ELSE
1739  ppatch%nBCellsVirt = 0
1740 
1741  nullify(ppatch%bvc)
1742  END IF ! pPatch%bcType
1743  END DO ! iPatch
1744 
1745 ! ******************************************************************************
1746 ! End
1747 ! ******************************************************************************
1748 
1749  IF ( global%myProcid == masterproc .AND. &
1750  global%verbLevel > verbose_low ) THEN
1751  WRITE(stdout,'(A,1X,A)') solver_name, &
1752  'Building patch-to-virtual-cell list done.'
1753  END IF ! global%verbLevel
1754 
1755  CALL deregisterfunction(global)
1756 
1757  END SUBROUTINE rflu_sype_buildp2vclistserial
1758 
1759 
1760 
1761 
1762 
1763 ! ******************************************************************************
1764 !
1765 ! Purpose: Build transforms between periodic patches.
1766 !
1767 ! Description: None.
1768 !
1769 ! Input:
1770 ! pRegion Pointer to region data
1771 !
1772 ! Output: None.
1773 !
1774 ! Notes: None.
1775 !
1776 ! ******************************************************************************
1777 
1778 SUBROUTINE rflu_sype_buildtransforms(pRegion)
1779 
1780  USE modtools
1781 
1782  IMPLICIT NONE
1783 
1784 ! ******************************************************************************
1785 ! Declarations and definitions
1786 ! ******************************************************************************
1787 
1788 ! ==============================================================================
1789 ! Arguments
1790 ! ==============================================================================
1791 
1792  TYPE(t_region), POINTER :: pregion
1793 
1794 ! ==============================================================================
1795 ! Locals
1796 ! ==============================================================================
1797 
1798  INTEGER :: ipatch,ipatchrelated
1799  REAL(RFREAL) :: ct,dotp,eqtol,ex,ey,ez,nx,nxr,ny,nyr,nz,nzr,st,theta
1800  TYPE(t_border), POINTER :: pborder
1801  TYPE(t_global), POINTER :: global
1802  TYPE(t_grid), POINTER :: pgrid
1803  TYPE(t_patch), POINTER :: ppatch,ppatchrelated
1804 
1805 ! ******************************************************************************
1806 ! Start
1807 ! ******************************************************************************
1808 
1809  global => pregion%global
1810 
1811  CALL registerfunction(global,'RFLU_SYPE_BuildTransforms',&
1812  'RFLU_ModSymmetryPeriodic.F90')
1813 
1814  IF ( global%myProcid == masterproc .AND. &
1815  global%verbLevel > verbose_none ) THEN
1816  WRITE(stdout,'(A,1X,A)') solver_name,'Building periodicity transforms...'
1817  END IF ! global%myProcid
1818 
1819 ! ******************************************************************************
1820 ! Set pointers and variables
1821 ! ******************************************************************************
1822 
1823  pgrid => pregion%grid
1824 
1825  eqtol = 1.0e-6_rfreal
1826 
1827 ! ******************************************************************************
1828 ! Loop over patches
1829 ! ******************************************************************************
1830 
1831  DO ipatch = 1,pgrid%nPatches
1832  ppatch => pregion%patches(ipatch)
1833 
1834 ! ==============================================================================
1835 ! Periodic patch, find related patch
1836 ! ==============================================================================
1837 
1838  IF ( ppatch%bcType == bc_periodic ) THEN
1839  ipatchrelated = ppatch%iPatchRelated
1840 
1841  ppatchrelated => pregion%patches(ipatchrelated)
1842 
1843 ! ------------------------------------------------------------------------------
1844 ! Check whether the two patches have the same flatness flag
1845 ! ------------------------------------------------------------------------------
1846 
1847  IF ( ppatch%flatFlag .EQV. ppatchrelated%flatFlag ) THEN
1848  nx = ppatch%pn(xcoord)
1849  ny = ppatch%pn(ycoord)
1850  nz = ppatch%pn(zcoord)
1851 
1852  nxr = ppatchrelated%pn(xcoord)
1853  nyr = ppatchrelated%pn(ycoord)
1854  nzr = ppatchrelated%pn(zcoord)
1855 
1856  dotp = nx*nxr + ny*nyr + nz*nzr
1857 
1858 ! ----- Linear periodicity -----------------------------------------------------
1859 
1860  IF ( (floatequal(dotp,-1.0_rfreal,eqtol) .EQV. .true.) ) THEN
1861  ppatch%tm(xcoord,xcoord) = 1.0_rfreal
1862  ppatch%tm(ycoord,xcoord) = 0.0_rfreal
1863  ppatch%tm(zcoord,xcoord) = 0.0_rfreal
1864  ppatch%tm(xyzmag,xcoord) = 0.0_rfreal
1865 
1866  ppatch%tm(xcoord,ycoord) = 0.0_rfreal
1867  ppatch%tm(ycoord,ycoord) = 1.0_rfreal
1868  ppatch%tm(zcoord,ycoord) = 0.0_rfreal
1869  ppatch%tm(xyzmag,ycoord) = 0.0_rfreal
1870 
1871  ppatch%tm(xcoord,zcoord) = 0.0_rfreal
1872  ppatch%tm(ycoord,zcoord) = 0.0_rfreal
1873  ppatch%tm(zcoord,zcoord) = 1.0_rfreal
1874  ppatch%tm(xyzmag,zcoord) = 0.0_rfreal
1875 
1876  ppatch%tm(xcoord,xyzmag) = ppatchrelated%pc(xcoord) - ppatch%pc(xcoord)
1877  ppatch%tm(ycoord,xyzmag) = ppatchrelated%pc(ycoord) - ppatch%pc(ycoord)
1878  ppatch%tm(zcoord,xyzmag) = ppatchrelated%pc(zcoord) - ppatch%pc(zcoord)
1879  ppatch%tm(xyzmag,xyzmag) = 1.0_rfreal
1880 
1881 ! ----- Rotational periodicity -------------------------------------------------
1882 
1883  ELSE
1884  theta = ppatch%angleRelated
1885 
1886  ct = cos(theta)
1887  st = sin(theta)
1888 
1889  SELECT CASE ( ppatch%axisRelated )
1890  CASE ( 1 )
1891  ex = 1.0_rfreal
1892  ey = 0.0_rfreal
1893  ez = 0.0_rfreal
1894  CASE ( 2 )
1895  ex = 0.0_rfreal
1896  ey = 1.0_rfreal
1897  ez = 0.0_rfreal
1898  CASE ( 3 )
1899  ex = 0.0_rfreal
1900  ey = 0.0_rfreal
1901  ez = 1.0_rfreal
1902  CASE default
1903  CALL errorstop(global,err_reached_default,__line__)
1904  END SELECT ! pPatch%axisRelated
1905 
1906  ppatch%tm(xcoord,xcoord) = ct + (1.0_rfreal-ct)*ex*ex
1907  ppatch%tm(ycoord,xcoord) = (1.0_rfreal-ct)*ex*ey - st*ez
1908  ppatch%tm(zcoord,xcoord) = (1.0_rfreal-ct)*ex*ez + st*ey
1909  ppatch%tm(xyzmag,xcoord) = 0.0_rfreal
1910 
1911  ppatch%tm(xcoord,ycoord) = (1.0_rfreal-ct)*ey*ex + st*ez
1912  ppatch%tm(ycoord,ycoord) = ct + (1.0_rfreal-ct)*ey*ey
1913  ppatch%tm(zcoord,ycoord) = (1.0_rfreal-ct)*ey*ez - st*ex
1914  ppatch%tm(xyzmag,ycoord) = 0.0_rfreal
1915 
1916  ppatch%tm(xcoord,zcoord) = (1.0_rfreal-ct)*ez*ex - st*ey
1917  ppatch%tm(ycoord,zcoord) = (1.0_rfreal-ct)*ez*ey + st*ex
1918  ppatch%tm(zcoord,zcoord) = ct + (1.0_rfreal-ct)*ez*ez
1919  ppatch%tm(xyzmag,zcoord) = 0.0_rfreal
1920 
1921  ppatch%tm(xcoord,xyzmag) = 0.0_rfreal
1922  ppatch%tm(ycoord,xyzmag) = 0.0_rfreal
1923  ppatch%tm(zcoord,xyzmag) = 0.0_rfreal
1924  ppatch%tm(xyzmag,xyzmag) = 1.0_rfreal
1925  END IF ! FloatEqual
1926 
1927 ! ------------------------------------------------------------------------------
1928 ! Invalid combination - related patches must have same flatness flag
1929 ! ------------------------------------------------------------------------------
1930 
1931  ELSE
1932  CALL errorstop(global,err_flatflag_inconsistent,__line__)
1933  END IF ! pPatch%flatFlag
1934  END IF ! pPatch%bcType
1935  END DO ! iPatch
1936 
1937 ! ******************************************************************************
1938 ! End
1939 ! ******************************************************************************
1940 
1941  IF ( global%myProcid == masterproc .AND. &
1942  global%verbLevel > verbose_none ) THEN
1943  WRITE(stdout,'(A,1X,A)') solver_name,'Building periodicity transforms done.'
1944  END IF ! global%myProcid
1945 
1946  CALL deregisterfunction(global)
1947 
1948 END SUBROUTINE rflu_sype_buildtransforms
1949 
1950 
1951 
1952 
1953 
1954 
1955 
1956 
1957 
1958 
1959 ! ******************************************************************************
1960 !
1961 ! Purpose: Build vertex maps which link local vertex lists on related periodic
1962 ! patches.
1963 !
1964 ! Description: Coordinates on current patch are transformed to related patch,
1965 ! then use Octree to accelerate search for matching vertex pairs.
1966 !
1967 ! Input:
1968 ! pRegion Pointer to region data
1969 !
1970 ! Output: None.
1971 !
1972 ! Notes: None.
1973 !
1974 ! ******************************************************************************
1975 
1976 SUBROUTINE rflu_sype_buildvertexmaps(pRegion)
1977 
1978  USE modtools
1979 
1980  USE rflu_modoctree
1981 
1982  IMPLICIT NONE
1983 
1984 ! ******************************************************************************
1985 ! Declarations and definitions
1986 ! ******************************************************************************
1987 
1988 ! ==============================================================================
1989 ! Arguments
1990 ! ==============================================================================
1991 
1992  TYPE(t_region), POINTER :: pregion
1993 
1994 ! ==============================================================================
1995 ! Locals
1996 ! ==============================================================================
1997 
1998  INTEGER :: cvsize,errorflag,icv,ipatch,ipatchrelated,ivg,ivg2,ivl,ivlmin,ivl2
1999  INTEGER, DIMENSION(:), ALLOCATABLE :: cv
2000  REAL(RFREAL) :: delfrac,dist,distmax,distmin,distminsave,eqtol,maxfrac, &
2001  xdel,xmax,xmin,ydel,ymax,ymin,zdel,zmax,zmin
2002  REAL(RFREAL) :: xyz(xcoord:xyzmag),xyzt(xcoord:xyzmag)
2003  REAL(RFREAL), DIMENSION(:,:), ALLOCATABLE :: xyzr
2004  TYPE(t_global), POINTER :: global
2005  TYPE(t_grid), POINTER :: pgrid
2006  TYPE(t_patch), POINTER :: ppatch,ppatchrelated
2007 
2008 ! ******************************************************************************
2009 ! Start
2010 ! ******************************************************************************
2011 
2012  global => pregion%global
2013 
2014  CALL registerfunction(global,'RFLU_SYPE_BuildVertexMaps',&
2015  'RFLU_ModSymmetryPeriodic.F90')
2016 
2017  IF ( global%myProcid == masterproc .AND. &
2018  global%verbLevel > verbose_none ) THEN
2019  WRITE(stdout,'(A,1X,A)') solver_name,'Building vertex maps...'
2020  END IF ! global%myProcid
2021 
2022 ! ******************************************************************************
2023 ! Set pointers and constants
2024 ! ******************************************************************************
2025 
2026  pgrid => pregion%grid
2027 
2028  cvsize = 10
2029 
2030  delfrac = 0.01_rfreal
2031  maxfrac = 0.1_rfreal
2032 
2033 ! ******************************************************************************
2034 ! Allocate temporary memory
2035 ! ******************************************************************************
2036 
2037  ALLOCATE(cv(cvsize),stat=errorflag)
2038  global%error = errorflag
2039  IF ( global%error /= err_none ) THEN
2040  CALL errorstop(global,err_allocate,__line__,'cv')
2041  END IF ! global%error
2042 
2043 ! ******************************************************************************
2044 ! Loop over patches
2045 ! ******************************************************************************
2046 
2047  DO ipatch = 1,pgrid%nPatches
2048  ppatch => pregion%patches(ipatch)
2049 
2050 ! ==============================================================================
2051 ! Set tolerance for vertex matching
2052 ! ==============================================================================
2053 
2054 
2055 ! ==============================================================================
2056 ! Periodic patch, find related patch
2057 ! ==============================================================================
2058 
2059  IF ( ppatch%bcType == bc_periodic ) THEN
2060  ipatchrelated = ppatch%iPatchRelated
2061  ppatchrelated => pregion%patches(ipatchrelated)
2062 
2063  IF ( global%myProcid == masterproc .AND. &
2064  global%verbLevel > verbose_none ) THEN
2065  WRITE(stdout,'(A,2X,2(1X,A,1X,I2),A)') solver_name,'Matching patches', &
2066  ipatch,'and',ipatchrelated,'...'
2067  END IF ! global%myProcid
2068 
2069 ! ------------------------------------------------------------------------------
2070 ! Index of related patch greater than that of current patch, so match
2071 ! vertices
2072 ! ------------------------------------------------------------------------------
2073 
2074  IF ( ipatchrelated > ipatch ) THEN
2075 
2076 ! ----- Initialize and write info ---------------------------------------------
2077 
2078  eqtol = 0.1_rfreal*sqrt(minval(ppatch%fn(xyzmag,1:ppatch%nBFaces)))
2079 
2080  distmax = -huge(1.0_rfreal)
2081  distminsave = huge(1.0_rfreal)
2082 
2083  IF ( global%myProcid == masterproc .AND. &
2084  global%verbLevel > verbose_none ) THEN
2085  WRITE(stdout,'(A,5X,A,24X,E13.6)') solver_name, &
2086  'Distance tolerance:',eqtol
2087  END IF ! global%myProcid
2088 
2089 ! ----- Allocate temporary memory ----------------------------------------------
2090 
2091  ALLOCATE(xyzr(xcoord:zcoord,ppatchrelated%nBVert),stat=errorflag)
2092  global%error = errorflag
2093  IF ( global%error /= err_none ) THEN
2094  CALL errorstop(global,err_allocate,__line__,'xyzr')
2095  END IF ! global%error
2096 
2097 ! ----- Set coordinates of vertices on related patch ---------------------------
2098 
2099  DO ivl = 1,ppatchrelated%nBVert
2100  ivg = ppatchrelated%bv(ivl)
2101 
2102  xyzr(xcoord,ivl) = pgrid%xyz(xcoord,ivg)
2103  xyzr(ycoord,ivl) = pgrid%xyz(ycoord,ivg)
2104  xyzr(zcoord,ivl) = pgrid%xyz(zcoord,ivg)
2105  END DO ! ivl
2106 
2107 ! ----- Build Octree based on vertices on related patch ------------------------
2108 
2109  xmin = minval(xyzr(xcoord,1:ppatchrelated%nBVert))
2110  xmax = maxval(xyzr(xcoord,1:ppatchrelated%nBVert))
2111  ymin = minval(xyzr(ycoord,1:ppatchrelated%nBVert))
2112  ymax = maxval(xyzr(ycoord,1:ppatchrelated%nBVert))
2113  zmin = minval(xyzr(zcoord,1:ppatchrelated%nBVert))
2114  zmax = maxval(xyzr(zcoord,1:ppatchrelated%nBVert))
2115 
2116  xdel = xmax - xmin
2117  ydel = ymax - ymin
2118  zdel = zmax - zmin
2119 
2120  xdel = max(xdel,min(maxfrac*ydel,maxfrac*zdel))
2121  ydel = max(ydel,min(maxfrac*xdel,maxfrac*zdel))
2122  zdel = max(zdel,min(maxfrac*xdel,maxfrac*ydel))
2123 
2124  xmin = xmin - delfrac*xdel
2125  xmax = xmax + delfrac*xdel
2126  ymin = ymin - delfrac*ydel
2127  ymax = ymax + delfrac*ydel
2128  zmin = zmin - delfrac*zdel
2129  zmax = zmax + delfrac*zdel
2130 
2131  CALL rflu_createoctree(global,ppatchrelated%nBVert)
2132  CALL rflu_buildoctree(xyzr(xcoord,1:ppatchrelated%nBVert), &
2133  xyzr(ycoord,1:ppatchrelated%nBVert), &
2134  xyzr(zcoord,1:ppatchrelated%nBVert), &
2136 
2137 ! ----- Loop over vertices on current patch ------------------------------------
2138 
2139  DO ivl = 1,ppatch%nBVert
2140  ivg = ppatch%bv(ivl)
2141 
2142 ! ------- Transform coordinates of current vertex
2143 
2144  xyz(xcoord) = pgrid%xyz(xcoord,ivg)
2145  xyz(ycoord) = pgrid%xyz(ycoord,ivg)
2146  xyz(zcoord) = pgrid%xyz(zcoord,ivg)
2147  xyz(xyzmag) = 1.0_rfreal
2148 
2149  xyzt = matmul(ppatch%tm,xyz)
2150 
2151 ! ------- Find closest <cvSize> vertices on related patch
2152 
2153  CALL rflu_queryoctree(xyzt(xcoord),xyzt(ycoord),xyzt(zcoord), &
2154  cvsize,cv)
2155 
2156 ! ------- Find closest vertex
2157 
2158  distmin = huge(1.0_rfreal)
2159 
2160  cvloop: DO icv = 1,cvsize
2161  ivl2 = cv(icv)
2162 
2163  ivg2 = ppatchrelated%bv(ivl2)
2164 
2165  dist = (xyzt(xcoord)-pgrid%xyz(xcoord,ivg2))**2 &
2166  + (xyzt(ycoord)-pgrid%xyz(ycoord,ivg2))**2 &
2167  + (xyzt(zcoord)-pgrid%xyz(zcoord,ivg2))**2
2168 
2169  IF ( dist > distmax ) THEN
2170  distmax = dist
2171  END IF ! dist
2172 
2173  IF ( dist <= distmin ) THEN
2174  distmin = dist
2175  ivlmin = ivl2
2176 
2177  IF ( distmin <= eqtol ) THEN
2178  EXIT cvloop
2179  END IF ! distMin
2180  END IF ! dist
2181  END DO cvloop
2182 
2183 ! ------- Store closest vertex in mapping array
2184 
2185  IF ( distmin <= eqtol ) THEN
2186  ppatch%bv2bv(ivl) = ivlmin
2187  ELSE
2188  CALL errorstop(global,err_vertex_match_failed,__line__)
2189  END IF ! distMin
2190 
2191  distminsave = min(distmin,distminsave)
2192  END DO ! ivl
2193 
2194 ! ----- Destroy Octree and deallocate temporary memory -------------------------
2195 
2196  CALL rflu_destroyoctree(global)
2197 
2198  DEALLOCATE(xyzr,stat=errorflag)
2199  global%error = errorflag
2200  IF ( global%error /= err_none ) THEN
2201  CALL errorstop(global,err_deallocate,__line__,'xyzr')
2202  END IF ! global%error
2203 
2204 ! ----- Write some info --------------------------------------------------------
2205 
2206  IF ( global%myProcid == masterproc .AND. &
2207  global%verbLevel > verbose_none ) THEN
2208  WRITE(stdout,'(A,5X,A,1X,E13.6)') solver_name, &
2209  'Maximum distance between matched vertices:',distmax
2210  WRITE(stdout,'(A,5X,A,1X,E13.6)') solver_name, &
2211  'Minimum distance between matched vertices:',distminsave
2212  END IF ! global%myProcid
2213 
2214 ! ------------------------------------------------------------------------------
2215 ! Index of related patch less than that of current patch, so copy previously
2216 ! built vertex matching array
2217 ! ------------------------------------------------------------------------------
2218 
2219  ELSE
2220  DO ivl = 1,ppatch%nBVert
2221  ivl2 = ppatchrelated%bv2bv(ivl)
2222 
2223  ppatch%bv2bv(ivl2) = ivl
2224  END DO ! ivl
2225  END IF ! iPatchRelated
2226 
2227  END IF ! pPatch%bcType
2228  END DO ! iPatch
2229 
2230 ! ******************************************************************************
2231 ! Deallocate temporary memory
2232 ! ******************************************************************************
2233 
2234  DEALLOCATE(cv,stat=errorflag)
2235  global%error = errorflag
2236  IF ( global%error /= err_none ) THEN
2237  CALL errorstop(global,err_deallocate,__line__,'cv')
2238  END IF ! global%error
2239 
2240 ! ******************************************************************************
2241 ! End
2242 ! ******************************************************************************
2243 
2244  IF ( global%myProcid == masterproc .AND. &
2245  global%verbLevel > verbose_none ) THEN
2246  WRITE(stdout,'(A,1X,A)') solver_name,'Building vertex maps done.'
2247  END IF ! global%myProcid
2248 
2249  CALL deregisterfunction(global)
2250 
2251 END SUBROUTINE rflu_sype_buildvertexmaps
2252 
2253 
2254 
2255 
2256 
2257 
2258 ! ******************************************************************************
2259 !
2260 ! Purpose: Create list of virtual cells adjacent to periodic and symmetry
2261 ! patches.
2262 !
2263 ! Description: None.
2264 !
2265 ! Input:
2266 ! pRegion Pointer to region
2267 !
2268 ! Output: None.
2269 !
2270 ! Notes: None.
2271 !
2272 ! ******************************************************************************
2273 
2274  SUBROUTINE rflu_sype_createp2vclist(pRegion)
2275 
2276  IMPLICIT NONE
2277 
2278 ! ******************************************************************************
2279 ! Declarations and definitions
2280 ! ******************************************************************************
2281 
2282 ! ==============================================================================
2283 ! Arguments
2284 ! ==============================================================================
2285 
2286  TYPE(t_region), POINTER :: pregion
2287 
2288 ! ==============================================================================
2289 ! Locals
2290 ! ==============================================================================
2291 
2292  INTEGER :: errorflag,ipatch
2293  TYPE(t_global), POINTER :: global
2294  TYPE(t_grid), POINTER :: pgrid
2295  TYPE(t_patch), POINTER :: ppatch
2296 
2297 ! ******************************************************************************
2298 ! Start
2299 ! ******************************************************************************
2300 
2301  global => pregion%global
2302 
2303  CALL registerfunction(global,'RFLU_SYPE_CreateP2VCList',&
2304  'RFLU_ModSymmetryPeriodic.F90')
2305 
2306  IF ( global%myProcid == masterproc .AND. &
2307  global%verbLevel > verbose_low ) THEN
2308  WRITE(stdout,'(A,1X,A)') solver_name, &
2309  'Creating patch-to-virtual-cell list...'
2310  WRITE(stdout,'(A,3X,A,1X,I5.5)') solver_name,'Global region:', &
2311  pregion%iRegionGlobal
2312  END IF ! global%verbLevel
2313 
2314 ! ******************************************************************************
2315 ! Set grid pointer
2316 ! ******************************************************************************
2317 
2318  pgrid => pregion%grid
2319 
2320 ! ******************************************************************************
2321 ! Allocate memory
2322 ! ******************************************************************************
2323 
2324  DO ipatch = 1,pgrid%nPatches
2325  ppatch => pregion%patches(ipatch)
2326 
2327  IF ( ppatch%nBCellsVirt > 0 ) THEN
2328  ALLOCATE(ppatch%bvc(ppatch%nBCellsVirt),stat=errorflag)
2329  global%error = errorflag
2330  IF ( global%error /= err_none ) THEN
2331  CALL errorstop(global,err_allocate,__line__,'pPatch%bvc')
2332  END IF ! global%error
2333  ELSE
2334  nullify(ppatch%bvc)
2335  END IF ! pPatch%nBCellsVirt
2336  END DO ! iPatch
2337 
2338 ! ******************************************************************************
2339 ! End
2340 ! ******************************************************************************
2341 
2342  IF ( global%myProcid == masterproc .AND. &
2343  global%verbLevel > verbose_low ) THEN
2344  WRITE(stdout,'(A,1X,A)') solver_name, &
2345  'Creating patch-to-virtual-cell list done.'
2346  END IF ! global%verbLevel
2347 
2348  CALL deregisterfunction(global)
2349 
2350  END SUBROUTINE rflu_sype_createp2vclist
2351 
2352 
2353 
2354 
2355 
2356 
2357 ! ******************************************************************************
2358 !
2359 ! Purpose: Create vertex maps.
2360 !
2361 ! Description: None.
2362 !
2363 ! Input:
2364 ! pRegion Pointer to region data
2365 !
2366 ! Output: None.
2367 !
2368 ! Notes: None.
2369 !
2370 ! ******************************************************************************
2371 
2372 SUBROUTINE rflu_sype_createvertexmaps(pRegion)
2373 
2374  IMPLICIT NONE
2375 
2376 ! ******************************************************************************
2377 ! Declarations and definitions
2378 ! ******************************************************************************
2379 
2380 ! ==============================================================================
2381 ! Arguments
2382 ! ==============================================================================
2383 
2384  TYPE(t_region), POINTER :: pregion
2385 
2386 ! ==============================================================================
2387 ! Locals
2388 ! ==============================================================================
2389 
2390  INTEGER :: errorflag,ipatch
2391  TYPE(t_global), POINTER :: global
2392  TYPE(t_grid), POINTER :: pgrid
2393  TYPE(t_patch), POINTER :: ppatch
2394 
2395 ! ******************************************************************************
2396 ! Start
2397 ! ******************************************************************************
2398 
2399  global => pregion%global
2400 
2401  CALL registerfunction(global,'RFLU_SYPE_CreateVertexMaps',&
2402  'RFLU_ModSymmetryPeriodic.F90')
2403 
2404  IF ( global%myProcid == masterproc .AND. &
2405  global%verbLevel > verbose_none ) THEN
2406  WRITE(stdout,'(A,1X,A)') solver_name,'Creating vertex maps...'
2407  END IF ! global%myProcid
2408 
2409  pgrid => pregion%grid
2410 
2411 ! ******************************************************************************
2412 ! Loop over patches and allocate memory
2413 ! ******************************************************************************
2414 
2415  DO ipatch = 1,pgrid%nPatches
2416  ppatch => pregion%patches(ipatch)
2417 
2418  IF ( ppatch%bcType == bc_periodic ) THEN
2419  ALLOCATE(ppatch%bv2bv(ppatch%nBVert),stat=errorflag)
2420  global%error = errorflag
2421  IF ( global%error /= err_none ) THEN
2422  CALL errorstop(global,err_allocate,__line__,'pPatch%bv2bv')
2423  END IF ! global%error
2424  ELSE
2425  nullify(ppatch%bv2bv)
2426  END IF ! pPatch%bcType
2427  END DO ! iPatch
2428 
2429 ! ******************************************************************************
2430 ! End
2431 ! ******************************************************************************
2432 
2433  IF ( global%myProcid == masterproc .AND. &
2434  global%verbLevel > verbose_none ) THEN
2435  WRITE(stdout,'(A,1X,A)') solver_name,'Creating vertex maps done.'
2436  END IF ! global%myProcid
2437 
2438  CALL deregisterfunction(global)
2439 
2440 END SUBROUTINE rflu_sype_createvertexmaps
2441 
2442 
2443 
2444 
2445 
2446 
2447 ! ******************************************************************************
2448 !
2449 ! Purpose: Destroy list of virtual cells adjacent to periodic and symmetry
2450 ! patches.
2451 !
2452 ! Description: None.
2453 !
2454 ! Input:
2455 ! pRegion Pointer to region
2456 !
2457 ! Output: None.
2458 !
2459 ! Notes: None.
2460 !
2461 ! ******************************************************************************
2462 
2463  SUBROUTINE rflu_sype_destroyp2vclist(pRegion)
2464 
2465  IMPLICIT NONE
2466 
2467 ! ******************************************************************************
2468 ! Declarations and definitions
2469 ! ******************************************************************************
2470 
2471 ! ==============================================================================
2472 ! Arguments
2473 ! ==============================================================================
2474 
2475  TYPE(t_region), POINTER :: pregion
2476 
2477 ! ==============================================================================
2478 ! Locals
2479 ! ==============================================================================
2480 
2481  INTEGER :: errorflag,ipatch
2482  TYPE(t_global), POINTER :: global
2483  TYPE(t_grid), POINTER :: pgrid
2484  TYPE(t_patch), POINTER :: ppatch
2485 
2486 ! ******************************************************************************
2487 ! Start
2488 ! ******************************************************************************
2489 
2490  global => pregion%global
2491 
2492  CALL registerfunction(global,'RFLU_SYPE_DestroyP2VCList',&
2493  'RFLU_ModSymmetryPeriodic.F90')
2494 
2495  IF ( global%myProcid == masterproc .AND. &
2496  global%verbLevel > verbose_low ) THEN
2497  WRITE(stdout,'(A,1X,A)') solver_name, &
2498  'Destroying patch-to-virtual-cell list...'
2499  END IF ! global%verbLevel
2500 
2501 ! ******************************************************************************
2502 ! Set grid pointer
2503 ! ******************************************************************************
2504 
2505  pgrid => pregion%grid
2506 
2507 ! ******************************************************************************
2508 ! Deallocate memory
2509 ! ******************************************************************************
2510 
2511  DO ipatch = 1,pgrid%nPatches
2512  ppatch => pregion%patches(ipatch)
2513 
2514  IF ( ppatch%nBCellsVirt > 0 ) THEN
2515  ppatch%nBCellsVirt = 0
2516 
2517  DEALLOCATE(ppatch%bvc,stat=errorflag)
2518  global%error = errorflag
2519  IF ( global%error /= err_none ) THEN
2520  CALL errorstop(global,err_deallocate,__line__,'pPatch%bvc')
2521  END IF ! global%error
2522  ELSE
2523  nullify(ppatch%bvc)
2524  END IF ! pPatch%nBCellsVirt
2525  END DO ! iPatch
2526 
2527 ! ******************************************************************************
2528 ! End
2529 ! ******************************************************************************
2530 
2531  IF ( global%myProcid == masterproc .AND. &
2532  global%verbLevel > verbose_low ) THEN
2533  WRITE(stdout,'(A,1X,A)') solver_name, &
2534  'Destroying patch-to-virtual-cell list done.'
2535  END IF ! global%verbLevel
2536 
2537  CALL deregisterfunction(global)
2538 
2539  END SUBROUTINE rflu_sype_destroyp2vclist
2540 
2541 
2542 
2543 
2544 
2545 
2546 
2547 
2548 ! ******************************************************************************
2549 !
2550 ! Purpose: For given virtual cell index on serial region, find corresponding
2551 ! actual cell index on serial region.
2552 !
2553 ! Description: None.
2554 !
2555 ! Input:
2556 ! pRegion Pointer to region
2557 ! icgs Cell index (must be virtual cell index)
2558 !
2559 ! Output:
2560 ! icgs2 Cell index (actual cell index)
2561 !
2562 ! Notes: None.
2563 !
2564 ! ******************************************************************************
2565 
2566 SUBROUTINE rflu_sype_getactualserialcell(pRegionSerial,icgs,icgs2)
2567 
2568  IMPLICIT NONE
2569 
2570 ! ******************************************************************************
2571 ! Declarations and definitions
2572 ! ******************************************************************************
2573 
2574 ! ==============================================================================
2575 ! Arguments
2576 ! ==============================================================================
2577 
2578  INTEGER, INTENT(IN) :: icgs
2579  INTEGER, INTENT(OUT) :: icgs2
2580  TYPE(t_region), POINTER :: pregionserial
2581 
2582 ! ==============================================================================
2583 ! Locals
2584 ! ==============================================================================
2585 
2586  LOGICAL :: foundflag
2587  INTEGER :: iborderserial,iborderserial2,iloc
2588  TYPE(t_border), POINTER :: pborderserial,pborderserial2
2589  TYPE(t_global), POINTER :: global
2590  TYPE(t_grid), POINTER :: pgridserial
2591 
2592 ! ******************************************************************************
2593 ! Start
2594 ! ******************************************************************************
2595 
2596  global => pregionserial%global
2597 
2598  CALL registerfunction(global,'RFLU_SYPE_GetActualSerialCell',&
2599  'RFLU_ModSymmetryPeriodic.F90')
2600 
2601 ! ******************************************************************************
2602 ! Set pointers
2603 ! ******************************************************************************
2604 
2605  pgridserial => pregionserial%grid
2606 
2607 ! ******************************************************************************
2608 ! Check region index and cell kind
2609 ! ******************************************************************************
2610 
2611  IF ( pregionserial%iRegionGlobal /= 0 ) THEN
2612  CALL errorstop(global,err_region_id_invalid,__line__)
2613  END IF ! pRegionSerial%iRegionGlobal
2614 
2615  IF ( (icgs <= pgridserial%nCells ) .OR. &
2616  (icgs > pgridserial%nCellsTot) ) THEN
2617  CALL errorstop(global,err_cell_kind_invalid,__line__)
2618  END IF ! icgs
2619 
2620 ! ******************************************************************************
2621 ! Loop over serial borders and search for cell
2622 ! ******************************************************************************
2623 
2624  borderloop2: DO iborderserial = 1,pgridserial%nBorders
2625  pborderserial => pgridserial%borders(iborderserial)
2626 
2627  CALL binarysearchinteger(pborderserial%icgRecv(1:pborderserial%nCellsRecv), &
2628  pborderserial%nCellsRecv,icgs,iloc)
2629 
2630  IF ( iloc /= element_not_found ) THEN
2631  foundflag = .true.
2632 
2633  iborderserial2 = pborderserial%iBorder
2634 
2635  IF ( iborderserial2 > pgridserial%nBorders ) THEN
2636  CALL errorstop(global,err_border_index_invalid,__line__)
2637  END IF ! iBorderSerial2
2638 
2639  pborderserial2 => pgridserial%borders(iborderserial2)
2640 
2641  icgs2 = pborderserial2%icgSend(iloc)
2642 
2643  IF ( icgs2 > pgridserial%nCells ) THEN
2644  CALL errorstop(global,err_cell_kind_invalid,__line__)
2645  END IF ! icgs2
2646 
2647  EXIT borderloop2
2648  END IF ! iLoc
2649  END DO borderloop2
2650 
2651  IF ( foundflag .EQV. .false. ) THEN
2652  CALL errorstop(global,err_cell_not_found,__line__)
2653  END IF ! foundFlag
2654 
2655 ! ******************************************************************************
2656 ! End
2657 ! ******************************************************************************
2658 
2659  CALL deregisterfunction(global)
2660 
2661 END SUBROUTINE rflu_sype_getactualserialcell
2662 
2663 
2664 
2665 
2666 
2667 
2668 
2669 
2670 
2671 ! ******************************************************************************
2672 !
2673 ! Purpose: Given serial vertex index, find index of related vertex in another
2674 ! region.
2675 !
2676 ! Description: None.
2677 !
2678 ! Input:
2679 ! pRegionSerial Pointer to serial region data
2680 ! pRegion Pointer to region data
2681 ! ivgs Serial vertex index
2682 !
2683 ! Output:
2684 ! ivg Related vertex index
2685 !
2686 ! Notes: None.
2687 !
2688 ! ******************************************************************************
2689 
2690 SUBROUTINE rflu_sype_getrelatedvertex(pRegionSerial,pRegion,ivgs,ivg)
2691 
2692  IMPLICIT NONE
2693 
2694 ! ******************************************************************************
2695 ! Declarations and definitions
2696 ! ******************************************************************************
2697 
2698 ! ==============================================================================
2699 ! Arguments
2700 ! ==============================================================================
2701 
2702  INTEGER, INTENT(IN) :: ivgs
2703  INTEGER, INTENT(OUT) :: ivg
2704  TYPE(t_region), POINTER :: pregion,pregionserial
2705 
2706 ! ==============================================================================
2707 ! Locals
2708 ! ==============================================================================
2709 
2710  LOGICAL :: exitborderloopflag,foundflagrecv,foundflagshared
2711  INTEGER :: errorflag,iborderserial,iborderserial2,iloc,ivgs2,ivl,ivls2
2712  INTEGER, DIMENSION(:), ALLOCATABLE :: ivgsharedsorted,sortkey
2713  TYPE(t_border), POINTER :: pborderserial,pborderserial2
2714  TYPE(t_global), POINTER :: global
2715  TYPE(t_grid), POINTER :: pgrid,pgridserial
2716 
2717 ! ******************************************************************************
2718 ! Start
2719 ! ******************************************************************************
2720 
2721  global => pregionserial%global
2722 
2723  CALL registerfunction(global,'RFLU_SYPE_GetRelatedVertex',&
2724  'RFLU_ModSymmetryPeriodic.F90')
2725 
2726 ! ******************************************************************************
2727 ! Set pointers and variables
2728 ! ******************************************************************************
2729 
2730  pgridserial => pregionserial%grid
2731  pgrid => pregion%grid
2732 
2733  foundflagrecv = .false.
2734  foundflagshared = .false.
2735 
2736 ! ******************************************************************************
2737 ! Check region index
2738 ! ******************************************************************************
2739 
2740  IF ( pregionserial%iRegionGlobal /= 0 ) THEN
2741  CALL errorstop(global,err_region_id_invalid,__line__)
2742  END IF ! pRegionSerial%iRegionGlobal
2743 
2744 ! ******************************************************************************
2745 ! Search for vertex in vertices to be shared
2746 ! ******************************************************************************
2747 
2748  exitborderloopflag = .false.
2749 
2750  borderloopshared: DO iborderserial = 1,pgridserial%nBorders
2751  pborderserial => pgridserial%borders(iborderserial)
2752 
2753  ALLOCATE(ivgsharedsorted(pborderserial%nVertShared),stat=errorflag)
2754  global%error = errorflag
2755  IF ( global%error /= err_none ) THEN
2756  CALL errorstop(global,err_allocate,__line__,'ivgSharedSorted')
2757  END IF ! global%error
2758 
2759  ALLOCATE(sortkey(pborderserial%nVertShared),stat=errorflag)
2760  global%error = errorflag
2761  IF ( global%error /= err_none ) THEN
2762  CALL errorstop(global,err_allocate,__line__,'sortKey')
2763  END IF ! global%error
2764 
2765  DO ivl = 1,pborderserial%nVertShared
2766  ivgsharedsorted(ivl) = pborderserial%ivgShared(ivl)
2767  sortkey(ivl) = ivl
2768  END DO ! ivl
2769 
2770  CALL quicksortintegerinteger(ivgsharedsorted,sortkey, &
2771  pborderserial%nVertShared)
2772 
2773  CALL binarysearchinteger(ivgsharedsorted,pborderserial%nVertShared,ivgs, &
2774  iloc)
2775 
2776  IF ( iloc /= element_not_found ) THEN
2777  iborderserial2 = pborderserial%iBorder
2778 
2779  pborderserial2 => pgridserial%borders(iborderserial2)
2780 
2781  ivls2 = sortkey(iloc)
2782  ivgs2 = pborderserial2%ivgShared(ivls2)
2783 
2784  CALL binarysearchinteger(pgrid%sv2pv(1:1,1:pgrid%nVertTot), &
2785  pgrid%nVertTot,ivgs2,iloc)
2786 
2787  IF ( iloc /= element_not_found ) THEN
2788  ivg = pgrid%sv2pv(2,iloc)
2789 
2790  IF ( ivg > pgrid%nVert ) THEN ! Not actual vertex
2791  CALL errorstop(global,err_vertex_kind_invalid,__line__)
2792  END IF ! ivg
2793 
2794  exitborderloopflag = .true.
2795  foundflagshared = .true.
2796  ELSE
2797  CALL errorstop(global,err_vertex_not_found,__line__)
2798  END IF ! iLoc
2799  END IF ! iLoc
2800 
2801  DEALLOCATE(ivgsharedsorted,stat=errorflag)
2802  global%error = errorflag
2803  IF ( global%error /= err_none ) THEN
2804  CALL errorstop(global,err_deallocate,__line__,'ivgSharedSorted')
2805  END IF ! global%error
2806 
2807  DEALLOCATE(sortkey,stat=errorflag)
2808  global%error = errorflag
2809  IF ( global%error /= err_none ) THEN
2810  CALL errorstop(global,err_deallocate,__line__,'sortKey')
2811  END IF ! global%error
2812 
2813  IF ( exitborderloopflag .EQV. .true. ) THEN
2814  EXIT borderloopshared
2815  END IF ! exitBorderLoopFlag
2816  END DO borderloopshared
2817 
2818 ! ******************************************************************************
2819 ! Search for vertex in vertices to be received
2820 ! ******************************************************************************
2821 
2822  IF ( foundflagshared .EQV. .false. ) THEN
2823  exitborderloopflag = .false.
2824 
2825  borderlooprecv: DO iborderserial = 1,pgridserial%nBorders
2826  pborderserial => pgridserial%borders(iborderserial)
2827 
2828  CALL binarysearchinteger(pborderserial%ivgRecv(1:pborderserial%nVertRecv), &
2829  pborderserial%nVertRecv,ivgs,iloc)
2830 
2831  IF ( iloc /= element_not_found ) THEN
2832  iborderserial2 = pborderserial%iBorder
2833 
2834  pborderserial2 => pgridserial%borders(iborderserial2)
2835 
2836  ivgs2 = pborderserial2%ivgSend(iloc)
2837 
2838  CALL binarysearchinteger(pgrid%sv2pv(1:1,1:pgrid%nVertTot), &
2839  pgrid%nVertTot,ivgs2,iloc)
2840 
2841  IF ( iloc /= element_not_found ) THEN
2842  ivg = pgrid%sv2pv(2,iloc)
2843 
2844  IF ( ivg > pgrid%nVert ) THEN ! Not actual vertex
2845  CALL errorstop(global,err_vertex_kind_invalid,__line__)
2846  END IF ! ivg
2847 
2848  exitborderloopflag = .true.
2849  foundflagrecv = .true.
2850  ELSE
2851  CALL errorstop(global,err_vertex_not_found,__line__)
2852  END IF ! iLoc
2853  END IF ! iLoc
2854 
2855  IF ( exitborderloopflag .EQV. .true. ) THEN
2856  EXIT borderlooprecv
2857  END IF ! exitBorderLoopFlag
2858  END DO borderlooprecv
2859  END IF ! foundFlagShared
2860 
2861 
2862  IF ( (foundflagshared .EQV. .false.) .AND. &
2863  (foundflagrecv .EQV. .false.) ) THEN
2864  CALL errorstop(global,err_vertex_not_found,__line__)
2865  END IF ! foundFlagShared
2866 
2867 ! ******************************************************************************
2868 ! End
2869 ! ******************************************************************************
2870 
2871  CALL deregisterfunction(global)
2872 
2873 END SUBROUTINE rflu_sype_getrelatedvertex
2874 
2875 
2876 
2877 
2878 
2879 
2880 
2881 
2882 
2883 
2884 ! ******************************************************************************
2885 !
2886 ! Purpose: Determine whether have symmetry and/or periodic patches.
2887 !
2888 ! Description: None.
2889 !
2890 ! Input:
2891 ! pRegion Pointer to region
2892 !
2893 ! Output: None.
2894 !
2895 ! Notes: None.
2896 !
2897 ! ******************************************************************************
2898 
2899 LOGICAL FUNCTION rflu_sype_havesypepatches(pRegion)
2900 
2901  IMPLICIT NONE
2902 
2903 ! ******************************************************************************
2904 ! Declarations and definitions
2905 ! ******************************************************************************
2906 
2907 ! ==============================================================================
2908 ! Arguments
2909 ! ==============================================================================
2910 
2911  TYPE(t_region), POINTER :: pregion
2912 
2913 ! ==============================================================================
2914 ! Locals
2915 ! ==============================================================================
2916 
2917  INTEGER :: errorflag,ipatch,nsypepatches
2918  TYPE(t_global), POINTER :: global
2919  TYPE(t_grid), POINTER :: pgrid
2920  TYPE(t_patch), POINTER :: ppatch
2921 
2922 ! ******************************************************************************
2923 ! Start
2924 ! ******************************************************************************
2925 
2926  global => pregion%global
2927 
2928  CALL registerfunction(global,'RFLU_SYPE_HaveSyPePatches',&
2929  'RFLU_ModSymmetryPeriodic.F90')
2930 
2931 ! ******************************************************************************
2932 ! Set pointers
2933 ! ******************************************************************************
2934 
2935  pgrid => pregion%grid
2936 
2937  nsypepatches = 0
2938 
2939 ! ******************************************************************************
2940 ! Determine whether have symmetry or periodic patches
2941 ! ******************************************************************************
2942 
2943  ipatchloop: DO ipatch = 1,pgrid%nPatches
2944  ppatch => pregion%patches(ipatch)
2945 
2946  IF ( ppatch%bcType == bc_symmetry .OR. ppatch%bcType == bc_periodic ) THEN
2947  nsypepatches = nsypepatches + 1
2948 
2949  EXIT ipatchloop
2950  END IF ! pPatch%bcType
2951  END DO ipatchloop
2952 
2953  IF ( nsypepatches > 0 ) THEN
2954  rflu_sype_havesypepatches = .true.
2955  ELSE
2956  rflu_sype_havesypepatches = .false.
2957  END IF ! nPatches
2958 
2959 ! ******************************************************************************
2960 ! End
2961 ! ******************************************************************************
2962 
2963  CALL deregisterfunction(global)
2964 
2965 END FUNCTION rflu_sype_havesypepatches
2966 
2967 
2968 
2969 
2970 
2971 
2972 
2973 
2974 
2975 ! ******************************************************************************
2976 !
2977 ! Purpose: Impose vertex maps.
2978 !
2979 ! Description: None.
2980 !
2981 ! Input:
2982 ! pRegion Pointer to region data
2983 ! pPatch Pointer to patch
2984 ! ivg Vertex index
2985 !
2986 ! Output:
2987 ! ivgm Vertex index
2988 !
2989 ! Notes: None.
2990 !
2991 ! ******************************************************************************
2992 
2993 SUBROUTINE rflu_sype_imposevertexmap(pRegion,pPatch,ivg,ivgm)
2994 
2995  IMPLICIT NONE
2996 
2997 ! ******************************************************************************
2998 ! Declarations and definitions
2999 ! ******************************************************************************
3000 
3001 ! ==============================================================================
3002 ! Arguments
3003 ! ==============================================================================
3004 
3005  INTEGER, INTENT(IN) :: ivg
3006  INTEGER, INTENT(OUT) :: ivgm
3007  TYPE(t_patch), POINTER :: ppatch
3008  TYPE(t_region), POINTER :: pregion
3009 
3010 ! ==============================================================================
3011 ! Locals
3012 ! ==============================================================================
3013 
3014  INTEGER :: errorflag,iloc,ipatchrelated,ivlm
3015  TYPE(t_global), POINTER :: global
3016  TYPE(t_patch), POINTER :: ppatchrelated
3017 
3018 ! ******************************************************************************
3019 ! Start
3020 ! ******************************************************************************
3021 
3022  global => pregion%global
3023 
3024  CALL registerfunction(global,'RFLU_SYPE_ImposeVertexMap',&
3025  'RFLU_ModSymmetryPeriodic.F90')
3026 
3027 ! ******************************************************************************
3028 ! Get related patch
3029 ! ******************************************************************************
3030 
3031  ipatchrelated = ppatch%iPatchRelated
3032  ppatchrelated => pregion%patches(ipatchrelated)
3033 
3034 ! ******************************************************************************
3035 ! Search for vertex and impose mapping
3036 ! ******************************************************************************
3037 
3038  CALL binarysearchinteger(ppatch%bv(1:ppatch%nBVert),ppatch%nBVert,ivg,iloc)
3039 
3040  IF ( iloc /= element_not_found ) THEN
3041  ivlm = ppatch%bv2bv(iloc)
3042  ivgm = ppatchrelated%bv(ivlm)
3043  ELSE
3044  CALL errorstop(global,err_vertex_not_found,__line__)
3045  END IF ! iLoc
3046 
3047 ! ******************************************************************************
3048 ! End
3049 ! ******************************************************************************
3050 
3051  CALL deregisterfunction(global)
3052 
3053 END SUBROUTINE rflu_sype_imposevertexmap
3054 
3055 
3056 
3057 
3058 
3059 
3060 ! ******************************************************************************
3061 !
3062 ! Purpose: Read and impose transforms between periodic patches.
3063 !
3064 ! Description: None.
3065 !
3066 ! Input:
3067 ! pRegion Pointer to region data
3068 !
3069 ! Output: None.
3070 !
3071 ! Notes: None.
3072 !
3073 ! ******************************************************************************
3074 
3075 SUBROUTINE rflu_sype_readtransforms(pRegion)
3076 
3078 
3079  IMPLICIT NONE
3080 
3081 ! ******************************************************************************
3082 ! Declarations and definitions
3083 ! ******************************************************************************
3084 
3085 ! ==============================================================================
3086 ! Arguments
3087 ! ==============================================================================
3088 
3089  TYPE(t_region), POINTER :: pregion
3090 
3091 ! ==============================================================================
3092 ! Locals
3093 ! ==============================================================================
3094 
3095  CHARACTER(CHRLEN) :: ifilename,sectionstring
3096  INTEGER :: errorflag,ifile,ipatch,ipatchglobal
3097  REAL(RFREAL) :: tm(xcoord:xyzmag,xcoord:xyzmag)
3098  TYPE(t_global), POINTER :: global
3099  TYPE(t_grid), POINTER :: pgrid
3100  TYPE(t_patch), POINTER :: ppatch
3101 
3102 ! ******************************************************************************
3103 ! Start
3104 ! ******************************************************************************
3105 
3106  global => pregion%global
3107 
3108  CALL registerfunction(global,'RFLU_SYPE_ReadTransforms',&
3109  'RFLU_ModSymmetryPeriodic.F90')
3110 
3111  IF ( global%myProcid == masterproc ) THEN
3112  IF ( global%verbLevel > verbose_none ) THEN
3113  WRITE(stdout,'(A,1X,A)') solver_name,'Reading periodicity transforms...'
3114 
3115  IF ( global%verbLevel > verbose_low ) THEN
3116  WRITE(stdout,'(A,3X,A,1X,I5.5)') solver_name,'Global region:', &
3117  pregion%iRegionGlobal
3118  END IF ! global%verbLevel
3119  END IF ! global%verbLevel
3120  END IF ! global%myProcid
3121 
3122 ! ******************************************************************************
3123 ! Set pointers and variables
3124 ! ******************************************************************************
3125 
3126  pgrid => pregion%grid
3127 
3128 ! ******************************************************************************
3129 ! Open file
3130 ! ******************************************************************************
3131 
3132  ifile = if_ptmatrix
3133 
3134  CALL buildfilenameplain(global,filedest_indir,'.ptm',ifilename)
3135 
3136  OPEN(ifile,file=ifilename,form="FORMATTED",status="OLD",iostat=errorflag)
3137  global%error = errorflag
3138  IF ( global%error /= err_none ) THEN
3139  CALL errorstop(global,err_file_open,__line__,trim(ifilename))
3140  END IF ! global%error
3141 
3142 ! ******************************************************************************
3143 ! Read header
3144 ! ******************************************************************************
3145 
3146  READ(ifile,'(A)') sectionstring
3147  IF ( trim(sectionstring) /= &
3148  '# ROCFLU periodic-transformation matrix file' ) THEN
3149  CALL errorstop(global,err_invalid_marker,__line__,sectionstring)
3150  END IF ! TRIM
3151 
3152 ! ******************************************************************************
3153 ! Read remainder of file
3154 ! ******************************************************************************
3155 
3156  emptyloop: DO
3157  READ(ifile,'(A)') sectionstring
3158 
3159  SELECT CASE ( trim(sectionstring) )
3160 
3161 ! ==============================================================================
3162 ! Patch section
3163 ! ==============================================================================
3164 
3165  CASE ( '# Patch' )
3166  READ(ifile,'(I3)') ipatchglobal
3167 
3168  IF ( global%myProcid == masterproc .AND. &
3169  global%verbLevel > verbose_low ) THEN
3170  WRITE(stdout,'(A,3X,A,1X,I3)') solver_name,'Global patch:', &
3171  ipatchglobal
3172  END IF ! global%myProcid
3173 
3174  READ(ifile,'(4(E23.16))') tm(xcoord,xcoord),tm(xcoord,ycoord), &
3175  tm(xcoord,zcoord),tm(xcoord,xyzmag)
3176  READ(ifile,'(4(E23.16))') tm(ycoord,xcoord),tm(ycoord,ycoord), &
3177  tm(ycoord,zcoord),tm(ycoord,xyzmag)
3178  READ(ifile,'(4(E23.16))') tm(zcoord,xcoord),tm(zcoord,ycoord), &
3179  tm(zcoord,zcoord),tm(zcoord,xyzmag)
3180  READ(ifile,'(4(E23.16))') tm(xyzmag,xcoord),tm(xyzmag,ycoord), &
3181  tm(xyzmag,zcoord),tm(xyzmag,xyzmag)
3182 
3183 ! ------------------------------------------------------------------------------
3184 ! Loop over local patches and store transformation matrix
3185 ! ------------------------------------------------------------------------------
3186 
3187  patchloop: DO ipatch = 1,pgrid%nPatches
3188  ppatch => pregion%patches(ipatch)
3189 
3190  IF ( ppatch%iPatchGlobal == ipatchglobal ) THEN
3191  IF ( global%myProcid == masterproc .AND. &
3192  global%verbLevel > verbose_low ) THEN
3193  WRITE(stdout,'(A,5X,A,1X,I3)') solver_name,'Local patch:', &
3194  ipatch,'.'
3195  END IF ! global%myProcid
3196 
3197  ppatch%transformFlag = .true. ! Set flag for later check
3198 
3199  ppatch%tm(xcoord,xcoord) = tm(xcoord,xcoord)
3200  ppatch%tm(xcoord,ycoord) = tm(xcoord,ycoord)
3201  ppatch%tm(xcoord,zcoord) = tm(xcoord,zcoord)
3202  ppatch%tm(xcoord,xyzmag) = tm(xcoord,xyzmag)
3203  ppatch%tm(ycoord,xcoord) = tm(ycoord,xcoord)
3204  ppatch%tm(ycoord,ycoord) = tm(ycoord,ycoord)
3205  ppatch%tm(ycoord,zcoord) = tm(ycoord,zcoord)
3206  ppatch%tm(ycoord,xyzmag) = tm(ycoord,xyzmag)
3207  ppatch%tm(zcoord,xcoord) = tm(zcoord,xcoord)
3208  ppatch%tm(zcoord,ycoord) = tm(zcoord,ycoord)
3209  ppatch%tm(zcoord,zcoord) = tm(zcoord,zcoord)
3210  ppatch%tm(zcoord,xyzmag) = tm(zcoord,xyzmag)
3211  ppatch%tm(xyzmag,xcoord) = tm(xyzmag,xcoord)
3212  ppatch%tm(xyzmag,ycoord) = tm(xyzmag,ycoord)
3213  ppatch%tm(xyzmag,zcoord) = tm(xyzmag,zcoord)
3214  ppatch%tm(xyzmag,xyzmag) = tm(xyzmag,xyzmag)
3215 
3216  EXIT patchloop
3217  END IF ! pPatch%iPatchGlobal
3218  END DO patchloop
3219 
3220 ! ==============================================================================
3221 ! Invalid section string
3222 ! ==============================================================================
3223 
3224  CASE ( '# End' )
3225  EXIT emptyloop
3226 
3227 ! ==============================================================================
3228 ! Default
3229 ! ==============================================================================
3230 
3231  CASE default
3232  IF ( global%verbLevel > verbose_low ) THEN
3233  WRITE(stdout,'(3X,A)') sectionstring
3234  END IF ! global%verbLevel
3235 
3236  CALL errorstop(global,err_invalid_marker,__line__,sectionstring)
3237  END SELECT ! TRIM
3238  END DO emptyloop
3239 
3240 ! ******************************************************************************
3241 ! Close file
3242 ! ******************************************************************************
3243 
3244  CLOSE(ifile,iostat=errorflag)
3245  global%error = errorflag
3246  IF ( global%error /= err_none ) THEN
3247  CALL errorstop(global,err_file_close,__line__,trim(ifilename))
3248  END IF ! global%error
3249 
3250 ! ******************************************************************************
3251 ! Check that periodic patches have transform
3252 ! ******************************************************************************
3253 
3254  DO ipatch = 1,pgrid%nPatches
3255  ppatch => pregion%patches(ipatch)
3256 
3257  IF ( ppatch%bcType == bc_periodic ) THEN
3258  IF ( ppatch%transformFlag .EQV. .false. ) THEN
3259 ! TEMPORARY
3260  WRITE(*,*) 'ERROR! Periodic patch without transform!'
3261  stop
3262 ! END TEMPORARY
3263  END IF ! pPatch%transformFlag
3264  END IF ! pPatch%bcType
3265  END DO ! iPatch
3266 
3267 ! ******************************************************************************
3268 ! End
3269 ! ******************************************************************************
3270 
3271  IF ( global%myProcid == masterproc .AND. &
3272  global%verbLevel > verbose_none ) THEN
3273  WRITE(stdout,'(A,1X,A)') solver_name,'Reading periodicity transforms done.'
3274  END IF ! global%myProcid
3275 
3276  CALL deregisterfunction(global)
3277 
3278 END SUBROUTINE rflu_sype_readtransforms
3279 
3280 
3281 
3282 
3283 
3284 
3285 
3286 
3287 ! ******************************************************************************
3288 !
3289 ! Purpose: Determine whether have symmetry and/or periodic patches by peeking
3290 ! into boundary-condition file without actually reading any information.
3291 !
3292 ! Description: None.
3293 !
3294 ! Input:
3295 ! global Pointer to global data
3296 !
3297 ! Output: None.
3298 !
3299 ! Notes:
3300 ! 1. Required because may need to know whether have patches before patch data
3301 ! structure is created, so cannot use RFLU_SYPE_HaveSyPePatches. This need
3302 ! occurs in rflupart when converting external grid representations to that
3303 ! used by Rocflu and need to estimate max number of faces, vertices, etc.
3304 ! Those max values will have to be larger than when do not have sype
3305 ! patches because virtual cells and vertices are going to be added.
3306 !
3307 ! ******************************************************************************
3308 
3309 SUBROUTINE rflu_sype_setsypepatchesflag(global)
3310 
3312 
3313  IMPLICIT NONE
3314 
3315 ! ******************************************************************************
3316 ! Declarations and definitions
3317 ! ******************************************************************************
3318 
3319 ! ==============================================================================
3320 ! Arguments
3321 ! ==============================================================================
3322 
3323  TYPE(t_global), POINTER :: global
3324 
3325 ! ==============================================================================
3326 ! Locals
3327 ! ==============================================================================
3328 
3329  CHARACTER(CHRLEN) :: ifilename
3330  CHARACTER(256) :: line
3331  INTEGER :: errorflag,ipatch,loopcounter
3332 
3333 ! ******************************************************************************
3334 ! Start
3335 ! ******************************************************************************
3336 
3337  CALL registerfunction(global,'RFLU_SYPE_SetSyPePatchesFlag',&
3338  'RFLU_ModSymmetryPeriodic.F90')
3339 
3340 ! ******************************************************************************
3341 ! Set variables
3342 ! ******************************************************************************
3343 
3344  loopcounter = 0
3345 
3346  global%syPePatchesFlag = .false.
3347 
3348 ! ******************************************************************************
3349 ! Open file
3350 ! ******************************************************************************
3351 
3352  CALL buildfilenameplain(global,filedest_indir,'.bc',ifilename)
3353 
3354  OPEN(if_input,file=ifilename,form='FORMATTED',status='OLD',iostat=errorflag)
3355  global%error = errorflag
3356  IF ( global%error /= err_none ) THEN
3357  CALL errorstop(global,err_file_open,__line__,'File: '//trim(ifilename))
3358  END IF ! global%error
3359 
3360 ! ******************************************************************************
3361 ! Look for
3362 ! ******************************************************************************
3363 
3364  keywordloop: DO
3365  READ(if_input,'(A256)',iostat=errorflag) line
3366 
3367  IF ( errorflag > 0 ) THEN ! Error occurred
3368  CALL errorstop(global,err_file_read,__line__,'File: '//trim(ifilename))
3369  ELSE IF ( errorflag < 0 ) THEN ! Encountered end of file
3370  EXIT keywordloop
3371  END IF ! errorFlag
3372 
3373  SELECT CASE( trim(line) )
3374  CASE ('# BC_PERIODIC')
3375  global%syPePatchesFlag = .true.
3376  EXIT keywordloop
3377  CASE ('# BC_SYMMETRY')
3378  global%syPePatchesFlag = .true.
3379  EXIT keywordloop
3380  CASE ('# END')
3381  EXIT keywordloop
3382  END SELECT ! TRIM(line)
3383 
3384  loopcounter = loopcounter + 1
3385 
3386  IF ( loopcounter >= limit_infinite_loop ) THEN ! Prevent infinite loop
3387  CALL errorstop(global,err_infinite_loop ,__line__)
3388  END IF ! loopCounter
3389  END DO keywordloop
3390 
3391 ! ******************************************************************************
3392 ! Close file
3393 ! ******************************************************************************
3394 
3395  CLOSE(if_input,iostat=errorflag)
3396  global%error = errorflag
3397  IF ( global%error /= err_none ) THEN
3398  CALL errorstop(global,err_file_close,__line__,'File: '//trim(ifilename))
3399  END IF ! global%error
3400 
3401 ! ******************************************************************************
3402 ! End
3403 ! ******************************************************************************
3404 
3405  CALL deregisterfunction(global)
3406 
3407 END SUBROUTINE rflu_sype_setsypepatchesflag
3408 
3409 
3410 
3411 
3412 
3413 
3414 
3415 ! ******************************************************************************
3416 !
3417 ! Purpose: Write transforms between periodic patches.
3418 !
3419 ! Description: None.
3420 !
3421 ! Input:
3422 ! pRegion Pointer to region data
3423 !
3424 ! Output: None.
3425 !
3426 ! Notes:
3427 ! 1. Should only be called for serial region.
3428 !
3429 ! ******************************************************************************
3430 
3431 SUBROUTINE rflu_sype_writetransforms(pRegion)
3432 
3434 
3435  IMPLICIT NONE
3436 
3437 ! ******************************************************************************
3438 ! Declarations and definitions
3439 ! ******************************************************************************
3440 
3441 ! ==============================================================================
3442 ! Arguments
3443 ! ==============================================================================
3444 
3445  TYPE(t_region), POINTER :: pregion
3446 
3447 ! ==============================================================================
3448 ! Locals
3449 ! ==============================================================================
3450 
3451  CHARACTER(CHRLEN) :: ifilename,sectionstring
3452  INTEGER :: errorflag,ifile,ipatch
3453  TYPE(t_global), POINTER :: global
3454  TYPE(t_grid), POINTER :: pgrid
3455  TYPE(t_patch), POINTER :: ppatch
3456 
3457 ! ******************************************************************************
3458 ! Start
3459 ! ******************************************************************************
3460 
3461  global => pregion%global
3462 
3463  CALL registerfunction(global,'RFLU_SYPE_WriteTransforms',&
3464  'RFLU_ModSymmetryPeriodic.F90')
3465 
3466  IF ( global%myProcid == masterproc ) THEN
3467  IF ( global%verbLevel > verbose_none ) THEN
3468  WRITE(stdout,'(A,1X,A)') solver_name,'Writing periodicity transforms...'
3469 
3470  IF ( global%verbLevel > verbose_low ) THEN
3471  WRITE(stdout,'(A,3X,A,1X,I5.5)') solver_name,'Global region:', &
3472  pregion%iRegionGlobal
3473  END IF ! global%verbLevel
3474  END IF ! global%verbLevel
3475  END IF ! global%myProcid
3476 
3477 ! ******************************************************************************
3478 ! Set pointers and variables
3479 ! ******************************************************************************
3480 
3481  pgrid => pregion%grid
3482 
3483 ! ******************************************************************************
3484 ! Open file
3485 ! ******************************************************************************
3486 
3487  ifile = if_ptmatrix
3488 
3489  CALL buildfilenameplain(global,filedest_indir,'.ptm',ifilename)
3490 
3491  OPEN(ifile,file=ifilename,form="FORMATTED",status="UNKNOWN",iostat=errorflag)
3492  global%error = errorflag
3493  IF ( global%error /= err_none ) THEN
3494  CALL errorstop(global,err_file_open,__line__,trim(ifilename))
3495  END IF ! global%error
3496 
3497 ! ******************************************************************************
3498 ! Write header
3499 ! ******************************************************************************
3500 
3501  sectionstring = '# ROCFLU periodic-transformation matrix file'
3502  WRITE(ifile,'(A)') trim(sectionstring)
3503 
3504 ! ******************************************************************************
3505 ! Loop over patches, write transformation matrix
3506 ! ******************************************************************************
3507 
3508  DO ipatch = 1,pgrid%nPatches
3509  ppatch => pregion%patches(ipatch)
3510 
3511  IF ( ppatch%bcType == bc_periodic ) THEN
3512  WRITE(ifile,'(A)') '# Patch'
3513 
3514  WRITE(ifile,'(I3)') ppatch%iPatchGlobal
3515 
3516  WRITE(ifile,'(4(E23.16))') ppatch%tm(xcoord,xcoord), &
3517  ppatch%tm(xcoord,ycoord), &
3518  ppatch%tm(xcoord,zcoord), &
3519  ppatch%tm(xcoord,xyzmag)
3520  WRITE(ifile,'(4(E23.16))') ppatch%tm(ycoord,xcoord), &
3521  ppatch%tm(ycoord,ycoord), &
3522  ppatch%tm(ycoord,zcoord), &
3523  ppatch%tm(ycoord,xyzmag)
3524  WRITE(ifile,'(4(E23.16))') ppatch%tm(zcoord,xcoord), &
3525  ppatch%tm(zcoord,ycoord), &
3526  ppatch%tm(zcoord,zcoord), &
3527  ppatch%tm(zcoord,xyzmag)
3528  WRITE(ifile,'(4(E23.16))') ppatch%tm(xyzmag,xcoord), &
3529  ppatch%tm(xyzmag,ycoord), &
3530  ppatch%tm(xyzmag,zcoord), &
3531  ppatch%tm(xyzmag,xyzmag)
3532  END IF ! pPatch%bcType
3533  END DO ! iPatch
3534 
3535 ! ******************************************************************************
3536 ! Write footer
3537 ! ******************************************************************************
3538 
3539  sectionstring = '# End'
3540  WRITE(ifile,'(A)') trim(sectionstring)
3541 
3542 ! ******************************************************************************
3543 ! Close file
3544 ! ******************************************************************************
3545 
3546  CLOSE(ifile,iostat=errorflag)
3547  global%error = errorflag
3548  IF ( global%error /= err_none ) THEN
3549  CALL errorstop(global,err_file_close,__line__,trim(ifilename))
3550  END IF ! global%error
3551 
3552 ! ******************************************************************************
3553 ! End
3554 ! ******************************************************************************
3555 
3556  IF ( global%myProcid == masterproc .AND. &
3557  global%verbLevel > verbose_none ) THEN
3558  WRITE(stdout,'(A,1X,A)') solver_name,'Writing periodicity transforms done.'
3559  END IF ! global%myProcid
3560 
3561  CALL deregisterfunction(global)
3562 
3563 END SUBROUTINE rflu_sype_writetransforms
3564 
3565 
3566 
3567 
3568 
3569 
3570 
3571 END MODULE rflu_modsymmetryperiodic
3572 
3573 ! ******************************************************************************
3574 !
3575 ! RCS Revision history:
3576 !
3577 ! $Log: RFLU_ModSymmetryPeriodic.F90,v $
3578 ! Revision 1.8 2008/12/06 08:44:24 mtcampbe
3579 ! Updated license.
3580 !
3581 ! Revision 1.7 2008/11/19 22:17:35 mtcampbe
3582 ! Added Illinois Open Source License/Copyright
3583 !
3584 ! Revision 1.6 2006/12/18 02:31:06 haselbac
3585 ! Bug fix: Inconsistent IFs in P2VC routines
3586 !
3587 ! Revision 1.5 2006/08/18 14:03:15 haselbac
3588 ! Added read/write routines for transform matrix
3589 !
3590 ! Revision 1.4 2006/04/17 19:57:05 haselbac
3591 ! Bug fix: Removed setting of pBorder in destroyign of patch-to-virtual cell list
3592 !
3593 ! Revision 1.3 2006/04/12 16:10:05 haselbac
3594 ! Increased nLayers for 2nd order virtual cells to 1
3595 !
3596 ! Revision 1.2 2006/04/07 15:19:20 haselbac
3597 ! Removed tabs
3598 !
3599 ! Revision 1.1 2006/03/25 21:38:54 haselbac
3600 ! Initial revision
3601 !
3602 ! ******************************************************************************
3603 
3604 
3605 
3606 
3607 
3608 
3609 
3610 
3611 
3612 
3613 
3614 
3615 
3616 
3617 
3618 
3619 
3620 
3621 
3622 
3623 
3624 
3625 
3626 
3627 
subroutine, public rflu_sype_setsypepatchesflag(global)
subroutine rflu_sype_imposevertexmap(pRegion, pPatch, ivg, ivgm)
subroutine, public rflu_sype_buildtransforms(pRegion)
subroutine, public rflu_createc2cstencilwrapper(pRegion)
double ymin() const
double xmax() const
subroutine, public rflu_queryoctree(XPT, YPT, ZPT, NUMP, NEIGHP)
subroutine rflu_sype_createp2vclist(pRegion)
CImg< T > & line(const unsigned int y0)
Get a line.
Definition: CImg.h:18421
subroutine, public rflu_buildfacevertlist(global, pGrid, fList, fListDim, vList, vListDimMax, vListDim, errorFlag)
subroutine, public rflu_sype_buildvertexmaps(pRegion)
Vector_n max(const Array_n_const &v1, const Array_n_const &v2)
Definition: Vector_n.h:354
double xmin() const
subroutine rflu_sype_augmentcelllists(pRegion, pPatch, vc, nCellsVirt)
subroutine registerfunction(global, funName, fileName)
Definition: ModError.F90:449
int status() const
Obtain the status of the attribute.
Definition: Attribute.h:240
double zmin() const
double sqrt(double d)
Definition: double.h:73
subroutine, public rflu_buildoctree(XI, YI, ZI, XLOW, XUPP, YLOW, YUPP, ZLOW, ZUPP)
subroutine quicksortinteger(a, n)
subroutine binarysearchinteger(a, n, v, i, j)
subroutine, public rflu_sype_destroyp2vclist(pRegion)
subroutine buildfilenameplain(global, dest, ext, fileName)
subroutine, public rflu_sype_getactualserialcell(pRegionSerial, icgs, icgs2)
NT & sin
subroutine rflu_sype_addvirtualcellsinv1(pRegion, pPatch, vc, nCellsVirtMax, nCellsVirt)
subroutine, public rflu_destroyc2cstencilwrapper(pRegion)
subroutine, public rflu_destroyoctree(global)
subroutine quicksortintegerinteger(a, b, n)
subroutine, private rflu_sype_addvirtualcellsinv2(pRegion, pPatch, vc, nCellsVirtMax, nCellsVirt)
blockLoc i
Definition: read.cpp:79
subroutine, public rflu_sype_getrelatedvertex(pRegionSerial, pRegion, ivgs, ivg)
subroutine, public rflu_sype_addvirtualcells(pRegion)
**********************************************************************Rocstar Simulation Suite Illinois Rocstar LLC All rights reserved ****Illinois Rocstar LLC IL **www illinoisrocstar com **sales illinoisrocstar com WITHOUT WARRANTY OF ANY **EXPRESS OR INCLUDING BUT NOT LIMITED TO THE WARRANTIES **OF FITNESS FOR A PARTICULAR PURPOSE AND **NONINFRINGEMENT IN NO EVENT SHALL THE CONTRIBUTORS OR **COPYRIGHT HOLDERS BE LIABLE FOR ANY DAMAGES OR OTHER WHETHER IN AN ACTION OF TORT OR **Arising OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE **USE OR OTHER DEALINGS WITH THE SOFTWARE **********************************************************************INTERFACE SUBROUTINE form
double zmax() const
subroutine reflectposition(nx, ny, nz, xc, yc, zc, xComp, yComp, zComp)
double ymax() const
subroutine, public rflu_sype_buildp2vclist(pRegion, pRegionSerial)
subroutine, public rflu_sype_createvertexmaps(pRegion)
subroutine, public rflu_setinfoc2cstencilwrapper(pRegion, orderNominal)
subroutine, public rflu_buildc2cstencilwrapper(pRegion, icgInput, constrInput)
Vector_n min(const Array_n_const &v1, const Array_n_const &v2)
Definition: Vector_n.h:346
subroutine rflu_sype_augmentvertexlists(pRegion, pGrid, pPatch, pBorder, pBorderRelated, x2v, icl, nXTot, ivgSendSorted, ivgRecvSorted, ivgSharedSorted)
j indices j
Definition: Indexing.h:6
subroutine swapintegers(a, b)
Definition: ModTools.F90:53
subroutine errorstop(global, errorCode, errorLine, addMessage)
Definition: ModError.F90:483
subroutine swaprfreals(a, b)
Definition: ModTools.F90:69
subroutine, public rflu_sype_writetransforms(pRegion)
subroutine, public rflu_buildvertcellnghblist(global, pGrid, vListOrig, vListOrigDim, nLayers, iReg, cList, cListDimMax, cListDim)
subroutine, public rflu_sype_readtransforms(pRegion)
subroutine, public rflu_sype_buildp2vclistserial(pRegion)
long double dist(long double *coord1, long double *coord2, int size)
static T_Key key
Definition: vinci_lass.c:76
subroutine deregisterfunction(global)
Definition: ModError.F90:469
LOGICAL function, public rflu_sype_havesypepatches(pRegion)
NT & cos
subroutine, public rflu_createoctree(global, nPoints)
LOGICAL function floatequal(a, b, tolIn)
Definition: ModTools.F90:99
IndexType nvert() const
Definition: Mesh.H:565