Rocstar  1.0
Rocstar multiphysics simulation application
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
RFLU_ModMergeRegions.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: Suite of routines to merge regions.
26 !
27 ! Description: None.
28 !
29 ! Notes: None.
30 !
31 ! ******************************************************************************
32 !
33 ! $Id: RFLU_ModMergeRegions.F90,v 1.10 2008/12/06 08:45:06 mtcampbe Exp $
34 !
35 ! Copyright: (c) 2004-2005 by the University of Illinois
36 !
37 ! ******************************************************************************
38 
40 
41  USE modglobal, ONLY: t_global
42  USE moddatatypes
43  USE modparameters
44  USE moderror
45  USE modbndpatch, ONLY: t_patch
46  USE moddatastruct, ONLY: t_level,t_region
47  USE modgrid, ONLY: t_grid
48  USE modmpi
49 
50  IMPLICIT NONE
51 
52  PRIVATE
53  PUBLIC :: rflu_merg_mergegrid, &
56 
57 ! ******************************************************************************
58 ! Declarations and definitions
59 ! ******************************************************************************
60 
61  CHARACTER(CHRLEN) :: &
62  RCSIdentString = '$RCSfile: RFLU_ModMergeRegions.F90,v $ $Revision: 1.10 $'
63 
64 
65 ! ******************************************************************************
66 ! Routines
67 ! ******************************************************************************
68 
69  CONTAINS
70 
71 
72 
73 
74 
75 
76 
77 
78 ! ******************************************************************************
79 !
80 ! Purpose: Merge grids.
81 !
82 ! Description: None.
83 !
84 ! Input:
85 ! pRegion Pointer to region
86 ! pRegionSerial Pointer to serial region
87 !
88 ! Output: None.
89 !
90 ! Notes: None.
91 !
92 ! ******************************************************************************
93 
94  SUBROUTINE rflu_merg_mergegrid(pRegion,pRegionSerial)
95 
96  IMPLICIT NONE
97 
98 ! ******************************************************************************
99 ! Declarations and definitions
100 ! ******************************************************************************
101 
102 ! ==============================================================================
103 ! Arguments
104 ! ==============================================================================
105 
106  TYPE(t_region), POINTER :: pregion,pregionserial
107 
108 ! ==============================================================================
109 ! Locals
110 ! ==============================================================================
111 
112  INTEGER :: icg,icg2,icl,icl2,ict2,ifl,ifl2,ipatch,ivg,ivl,ivg2,offs
113  TYPE(t_grid), POINTER :: pgrid,pgridserial
114  TYPE(t_global), POINTER :: global
115  TYPE(t_patch), POINTER :: ppatch,ppatchserial
116 
117 ! ******************************************************************************
118 ! Start
119 ! ******************************************************************************
120 
121  global => pregion%global
122 
123  CALL registerfunction(global,'RFLU_MERG_MergeGrid',&
124  'RFLU_ModMergeRegions.F90')
125 
126  IF ( global%verbLevel > verbose_none ) THEN
127  WRITE(stdout,'(A,1X,A)') solver_name,'Merging grid...'
128  END IF ! global%verbLevel
129 
130  IF ( global%verbLevel > verbose_low ) THEN
131  WRITE(stdout,'(A,3X,A,1X,I5.5)') solver_name,'Global region:', &
132  pregion%iRegionGlobal
133  END IF ! global%verbLevel
134 
135 ! ******************************************************************************
136 ! Set pointers
137 ! ******************************************************************************
138 
139  pgrid => pregion%grid
140  pgridserial => pregionserial%grid
141 
142 ! ******************************************************************************
143 ! Vertex coordinates
144 ! ******************************************************************************
145 
146  IF ( global%verbLevel > verbose_low ) THEN
147  WRITE(stdout,'(A,3X,A)') solver_name,'Merging vertex coordinates...'
148  END IF ! global%verbLevel
149 
150  DO ivg = 1,pgrid%nVertTot
151  ivg2 = pgrid%pv2sv(ivg)
152 
153  pgridserial%xyz(xcoord,ivg2) = pgrid%xyz(xcoord,ivg)
154  pgridserial%xyz(ycoord,ivg2) = pgrid%xyz(ycoord,ivg)
155  pgridserial%xyz(zcoord,ivg2) = pgrid%xyz(zcoord,ivg)
156  END DO ! ivg
157 
158 ! ******************************************************************************
159 ! Cell connectivity
160 ! ******************************************************************************
161 
162  IF ( global%verbLevel > verbose_low ) THEN
163  WRITE(stdout,'(A,3X,A)') solver_name,'Merging cell connectivity...'
164  END IF ! global%verbLevel
165 
166 ! ==============================================================================
167 ! Tetrahedra
168 ! ==============================================================================
169 
170  IF ( (global%verbLevel > verbose_low) .AND. (pgrid%nTetsTot > 0) ) THEN
171  WRITE(stdout,'(A,5X,A)') solver_name,'Tetrahedra...'
172  END IF ! global%verbLevel
173 
174  DO icl = 1,pgrid%nTetsTot
175  icg = pgrid%tet2CellGlob(icl)
176  icg2 = pgrid%pc2sc(icg)
177 
178  ict2 = pgridserial%cellGlob2Loc(1,icg2)
179  icl2 = pgridserial%cellGlob2Loc(2,icg2)
180 
181  IF ( ict2 /= cell_type_tet ) THEN
182  CALL errorstop(global,err_cell_type_invalid,__line__)
183  END IF ! ict2
184 
185  DO ivl = 1,4
186  ivg = pgrid%tet2v(ivl,icl)
187  ivg2 = pgrid%pv2sv(ivg)
188 
189  pgridserial%tet2v(ivl,icl2) = ivg2
190  END DO ! ivl
191  END DO ! icl
192 
193 ! ==============================================================================
194 ! Hexahedra
195 ! ==============================================================================
196 
197  IF ( (global%verbLevel > verbose_low) .AND. (pgrid%nHexsTot > 0) ) THEN
198  WRITE(stdout,'(A,5X,A)') solver_name,'Hexahedra...'
199  END IF ! global%verbLevel
200 
201  DO icl = 1,pgrid%nHexsTot
202  icg = pgrid%hex2CellGlob(icl)
203  icg2 = pgrid%pc2sc(icg)
204 
205  ict2 = pgridserial%cellGlob2Loc(1,icg2)
206  icl2 = pgridserial%cellGlob2Loc(2,icg2)
207 
208  IF ( ict2 /= cell_type_hex ) THEN
209  CALL errorstop(global,err_cell_type_invalid,__line__)
210  END IF ! ict2
211 
212  DO ivl = 1,8
213  ivg = pgrid%hex2v(ivl,icl)
214  ivg2 = pgrid%pv2sv(ivg)
215 
216  pgridserial%hex2v(ivl,icl2) = ivg2
217  END DO ! ivl
218  END DO ! icl
219 
220 ! ==============================================================================
221 ! Prisms
222 ! ==============================================================================
223 
224  IF ( (global%verbLevel > verbose_low) .AND. (pgrid%nPrisTot > 0) ) THEN
225  WRITE(stdout,'(A,5X,A)') solver_name,'Prisms...'
226  END IF ! global%verbLevel
227 
228  DO icl = 1,pgrid%nPrisTot
229  icg = pgrid%pri2CellGlob(icl)
230  icg2 = pgrid%pc2sc(icg)
231 
232  ict2 = pgridserial%cellGlob2Loc(1,icg2)
233  icl2 = pgridserial%cellGlob2Loc(2,icg2)
234 
235  IF ( ict2 /= cell_type_pri ) THEN
236  CALL errorstop(global,err_cell_type_invalid,__line__)
237  END IF ! ict2
238 
239  DO ivl = 1,6
240  ivg = pgrid%pri2v(ivl,icl)
241  ivg2 = pgrid%pv2sv(ivg)
242 
243  pgridserial%pri2v(ivl,icl2) = ivg2
244  END DO ! ivl
245  END DO ! icl
246 
247 ! ==============================================================================
248 ! Pyramids
249 ! ==============================================================================
250 
251  IF ( (global%verbLevel > verbose_low) .AND. (pgrid%nPyrsTot > 0) ) THEN
252  WRITE(stdout,'(A,5X,A)') solver_name,'Pyramids...'
253  END IF ! global%verbLevel
254 
255  DO icl = 1,pgrid%nPyrsTot
256  icg = pgrid%pyr2CellGlob(icl)
257  icg2 = pgrid%pc2sc(icg)
258 
259  ict2 = pgridserial%cellGlob2Loc(1,icg2)
260  icl2 = pgridserial%cellGlob2Loc(2,icg2)
261 
262  IF ( ict2 /= cell_type_pyr ) THEN
263  CALL errorstop(global,err_cell_type_invalid,__line__)
264  END IF ! ict2
265 
266  DO ivl = 1,5
267  ivg = pgrid%pyr2v(ivl,icl)
268  ivg2 = pgrid%pv2sv(ivg)
269 
270  pgridserial%pyr2v(ivl,icl2) = ivg2
271  END DO ! ivl
272  END DO ! icl
273 
274 ! ******************************************************************************
275 ! Patch connectivity
276 ! ******************************************************************************
277 
278  IF ( global%myProcid == masterproc .AND. &
279  global%verbLevel > verbose_low ) THEN
280  WRITE(stdout,'(A,3X,A)') solver_name,'Merging patch connectivity...'
281  END IF ! global%verbLevel
282 
283 ! ==============================================================================
284 ! Loop over patches
285 ! ==============================================================================
286 
287  DO ipatch = 1,pgrid%nPatches
288  ppatch => pregion%patches(ipatch)
289  ppatchserial => pregionserial%patches(ppatch%iPatchGlobal)
290 
291  IF ( global%verbLevel > verbose_low ) THEN
292  WRITE(stdout,'(A,5X,A,2(1X,I2))') solver_name,'Patch:',ipatch, &
293  ppatch%iPatchGlobal
294  END IF ! global%verbLevel
295 
296 ! ------------------------------------------------------------------------------
297 ! Triangles
298 ! ------------------------------------------------------------------------------
299 
300  IF ( (global%verbLevel > verbose_low) .AND. (ppatch%nBTrisTot > 0) ) THEN
301  WRITE(stdout,'(A,7X,A)') solver_name,'Triangles...'
302  END IF ! global%verbLevel
303 
304  offs = pgrid%pbf2sbfCSRInfo(ipatch) - 1
305 
306  DO ifl = 1,ppatch%nBTrisTot
307  ifl2 = pgrid%pbf2sbfCSR(offs+ifl)
308 
309  DO ivl = 1,3
310  ivg = ppatch%bTri2v(ivl,ifl)
311  ivg2 = pgrid%pv2sv(ivg)
312 
313  ppatchserial%bTri2v(ivl,ifl2) = ivg2
314  END DO ! ivl
315  END DO ! ifl
316 
317 ! ------------------------------------------------------------------------------
318 ! Quadrilaterals
319 ! ------------------------------------------------------------------------------
320 
321  IF ( (global%verbLevel > verbose_low) .AND. (ppatch%nBQuadsTot > 0) ) THEN
322  WRITE(stdout,'(A,7X,A)') solver_name,'Quadrilaterals...'
323  END IF ! global%verbLevel
324 
325  offs = offs + ppatch%nBTrisTot
326 
327  DO ifl = 1,ppatch%nBQuadsTot
328  ifl2 = pgrid%pbf2sbfCSR(offs+ifl) - ppatchserial%nBTrisTot
329 
330  DO ivl = 1,4
331  ivg = ppatch%bQuad2v(ivl,ifl)
332  ivg2 = pgrid%pv2sv(ivg)
333 
334  ppatchserial%bQuad2v(ivl,ifl2) = ivg2
335  END DO ! ivl
336  END DO ! ifl
337  END DO ! iPatch
338 
339 ! ******************************************************************************
340 ! End
341 ! ******************************************************************************
342 
343  IF ( global%verbLevel > verbose_none ) THEN
344  WRITE(stdout,'(A,1X,A)') solver_name,'Merging grid done.'
345  END IF ! global%verbLevel
346 
347  CALL deregisterfunction(global)
348 
349  END SUBROUTINE rflu_merg_mergegrid
350 
351 
352 
353 
354 
355 
356 
357 
358 ! ******************************************************************************
359 !
360 ! Purpose: Merge patch coefficients.
361 !
362 ! Description: None.
363 !
364 ! Input:
365 ! pRegion Pointer to region
366 ! pRegionSerial Pointer to serial region
367 !
368 ! Output: None.
369 !
370 ! Notes: None.
371 !
372 ! ******************************************************************************
373 
374  SUBROUTINE rflu_merg_mergepatchcoeffs(pRegion,pRegionSerial)
375 
376  IMPLICIT NONE
377 
378 ! ******************************************************************************
379 ! Declarations and definitions
380 ! ******************************************************************************
381 
382 ! ==============================================================================
383 ! Arguments
384 ! ==============================================================================
385 
386  TYPE(t_region), POINTER :: pregion,pregionserial
387 
388 ! ==============================================================================
389 ! Locals
390 ! ==============================================================================
391 
392  INTEGER :: ifl,ifl2,ipatch,offs
393  TYPE(t_global), POINTER :: global
394  TYPE(t_grid), POINTER :: pgrid,pgridserial
395  TYPE(t_patch), POINTER :: ppatch,ppatchserial
396 
397 ! ******************************************************************************
398 ! Start
399 ! ******************************************************************************
400 
401  global => pregion%global
402 
403  CALL registerfunction(global,'RFLU_MERG_MergePatchCoeffs',&
404  'RFLU_ModMergeRegions.F90')
405 
406  IF ( global%verbLevel > verbose_none ) THEN
407  WRITE(stdout,'(A,1X,A)') solver_name,'Merging patch coefficients...'
408  END IF ! global%verbLevel
409 
410  IF ( global%verbLevel > verbose_low ) THEN
411  WRITE(stdout,'(A,3X,A,1X,I5.5)') solver_name,'Global region:', &
412  pregion%iRegionGlobal
413  END IF ! global%verbLevel
414 
415 ! ******************************************************************************
416 ! Set pointers
417 ! ******************************************************************************
418 
419  pgrid => pregion%grid
420  pgridserial => pregionserial%grid
421 
422 ! ******************************************************************************
423 ! Merge patch coefficients
424 ! ******************************************************************************
425 
426  DO ipatch = 1,pgrid%nPatches
427  ppatch => pregion%patches(ipatch)
428 
429  ppatchserial => pregionserial%patches(ppatch%iPatchGlobal)
430 
431  offs = pgrid%pbf2sbfCSRInfo(ipatch) - 1
432 
433  DO ifl = 1,ppatch%nBTris
434  ifl2 = pgrid%pbf2sbfCSR(offs+ifl)
435 
436  ppatchserial%cp( ifl2) = ppatch%cp( ifl)
437  ppatchserial%cf(xcoord,ifl2) = ppatch%cf(xcoord,ifl)
438  ppatchserial%cf(ycoord,ifl2) = ppatch%cf(ycoord,ifl)
439  ppatchserial%cf(zcoord,ifl2) = ppatch%cf(zcoord,ifl)
440  ppatchserial%ch( ifl2) = ppatch%ch( ifl)
441  END DO ! ifl
442 
443  offs = offs + ppatch%nBTrisTot
444 
445  DO ifl = 1,ppatch%nBQuads
446  ifl2 = pgrid%pbf2sbfCSR(offs+ifl)
447 
448  ppatchserial%cp( ifl2) = ppatch%cp( ifl)
449  ppatchserial%cf(xcoord,ifl2) = ppatch%cf(xcoord,ifl)
450  ppatchserial%cf(ycoord,ifl2) = ppatch%cf(ycoord,ifl)
451  ppatchserial%cf(zcoord,ifl2) = ppatch%cf(zcoord,ifl)
452  ppatchserial%ch( ifl2) = ppatch%ch( ifl)
453  END DO ! ifl
454  END DO ! iPatch
455 
456 ! ******************************************************************************
457 ! End
458 ! ******************************************************************************
459 
460  IF ( global%verbLevel > verbose_none ) THEN
461  WRITE(stdout,'(A,1X,A)') solver_name,'Merging patch coefficients done.'
462  END IF ! global%verbLevel
463 
464  CALL deregisterfunction(global)
465 
466  END SUBROUTINE rflu_merg_mergepatchcoeffs
467 
468 
469 
470 
471 
472 
473 
474 
475 
476 ! ******************************************************************************
477 !
478 ! Purpose: Merge solution
479 !
480 ! Description: None.
481 !
482 ! Input:
483 ! pRegion Pointer to region
484 ! pRegionSerial Pointer to serial region
485 !
486 ! Output: None.
487 !
488 ! Notes: None.
489 !
490 ! ******************************************************************************
491 
492  SUBROUTINE rflu_merg_mergesolwrapper(pRegion,pRegionSerial)
493 
496 
497 #ifdef PLAG
498  USE modpartlag, ONLY: t_plag
500 #endif
501 
502  IMPLICIT NONE
503 
504 ! ******************************************************************************
505 ! Declarations and definitions
506 ! ******************************************************************************
507 
508 ! ==============================================================================
509 ! Arguments
510 ! ==============================================================================
511 
512  TYPE(t_region), POINTER :: pregion,pregionserial
513 
514 ! ==============================================================================
515 ! Locals
516 ! ==============================================================================
517 
518  INTEGER :: icg,icg2
519  REAL(RFREAL), DIMENSION(:,:), POINTER :: pcv,pcvserial
520  TYPE(t_grid), POINTER :: pgrid
521  TYPE(t_global), POINTER :: global
522 #ifdef PLAG
523  TYPE(t_plag), POINTER :: pplag,pplagserial
524 #endif
525 
526 ! ******************************************************************************
527 ! Start
528 ! ******************************************************************************
529 
530  global => pregion%global
531 
532  CALL registerfunction(global,'RFLU_MERG_MergeSolWrapper',&
533  'RFLU_ModMergeRegions.F90')
534 
535  IF ( global%verbLevel > verbose_none ) THEN
536  WRITE(stdout,'(A,1X,A)') solver_name,'Merging solution...'
537  END IF ! global%verbLevel
538 
539  IF ( global%verbLevel > verbose_low ) THEN
540  WRITE(stdout,'(A,3X,A,1X,I5.5)') solver_name,'Global region:', &
541  pregion%iRegionGlobal
542  END IF ! global%verbLevel
543 
544 ! ******************************************************************************
545 ! Set pointers
546 ! ******************************************************************************
547 
548  pgrid => pregion%grid
549 
550 ! ******************************************************************************
551 ! Merge solutions
552 ! ******************************************************************************
553 
554 ! ==============================================================================
555 ! Mixture
556 ! ==============================================================================
557 
558  pregionserial%mixt%cvState = cv_mixt_state_cons
559 
560  CALL rflu_copy_celldatap2s_r2d(global,pgrid,pregion%mixt%cv, &
561  pregionserial%mixt%cv)
562 
563 ! ==============================================================================
564 ! Physical modules
565 ! ==============================================================================
566 
567 #ifdef SPEC
568 ! ------------------------------------------------------------------------------
569 ! Species
570 ! ------------------------------------------------------------------------------
571 
572  IF ( global%specUsed .EQV. .true. ) THEN
573  pregionserial%spec%cvState = cv_mixt_state_cons
574 
575  CALL rflu_copy_celldatap2s_r2d(global,pgrid,pregion%spec%cv, &
576  pregionserial%spec%cv)
577 
578  IF ( pregion%specInput%nSpeciesEE > 0 ) THEN
579  CALL rflu_copy_celldatap2s_r3d(global,pgrid,pregion%spec%eev, &
580  pregionserial%spec%eev)
581  END IF ! pRegion%specInput%nSpeciesEE
582  END IF ! global%specUsed
583 #endif
584 
585 #ifdef PLAG
586 ! ------------------------------------------------------------------------------
587 ! Particles
588 ! ------------------------------------------------------------------------------
589 
590  IF ( global%plagUsed .EQV. .true. ) THEN
591  pplag => pregion%plag
592  pplagserial => pregionserial%plag
593 
594  CALL plag_dstr_mergeparticlewrapper(pregion,pplag,pplagserial)
595  END IF ! global%plagUsed
596 #endif
597 
598 ! ******************************************************************************
599 ! End
600 ! ******************************************************************************
601 
602  IF ( global%verbLevel > verbose_none ) THEN
603  WRITE(stdout,'(A,1X,A)') solver_name,'Merging solution done.'
604  END IF ! global%verbLevel
605 
606  CALL deregisterfunction(global)
607 
608  END SUBROUTINE rflu_merg_mergesolwrapper
609 
610 
611 
612 
613 
614 
615 
616 
617 
618 
619 ! ******************************************************************************
620 ! End
621 ! ******************************************************************************
622 
623 END MODULE rflu_modmergeregions
624 
625 
626 ! ******************************************************************************
627 !
628 ! RCS Revision history:
629 !
630 ! $Log: RFLU_ModMergeRegions.F90,v $
631 ! Revision 1.10 2008/12/06 08:45:06 mtcampbe
632 ! Updated license.
633 !
634 ! Revision 1.9 2008/11/19 22:18:16 mtcampbe
635 ! Added Illinois Open Source License/Copyright
636 !
637 ! Revision 1.8 2007/03/12 23:33:02 haselbac
638 ! Prepared code for merging of particles
639 !
640 ! Revision 1.7 2006/12/18 02:32:34 haselbac
641 ! Bug fixes: Needed nXTots to get merging to work with SYPE
642 !
643 ! Revision 1.6 2006/04/07 15:19:26 haselbac
644 ! Removed tabs
645 !
646 ! Revision 1.5 2005/11/27 02:00:57 haselbac
647 ! Bug fix: Missing cvState, added merging for EEv
648 !
649 ! Revision 1.4 2005/08/19 02:36:02 haselbac
650 ! Adapted to changes in RFLU_ModCopyData
651 !
652 ! Revision 1.3 2005/04/15 15:08:33 haselbac
653 ! Added routine to merge patch coeffs, general clean-up
654 !
655 ! Revision 1.2 2005/01/17 19:52:07 haselbac
656 ! Removed debug statements
657 !
658 ! Revision 1.1 2004/12/29 20:58:54 haselbac
659 ! Initial revision
660 !
661 ! ******************************************************************************
662 
663 
664 
665 
666 
667 
668 
669 
670 
subroutine, public rflu_copy_celldatap2s_r3d(global, pGrid, var, varSerial)
subroutine registerfunction(global, funName, fileName)
Definition: ModError.F90:449
subroutine, public rflu_copy_celldatap2s_r2d(global, pGrid, var, varSerial)
subroutine, public rflu_merg_mergegrid(pRegion, pRegionSerial)
subroutine, public plag_dstr_mergeparticlewrapper(pRegion, pPlag, pPlag2)
subroutine, public rflu_merg_mergesolwrapper(pRegion, pRegionSerial)
subroutine, public rflu_merg_mergepatchcoeffs(pRegion, pRegionSerial)
subroutine errorstop(global, errorCode, errorLine, addMessage)
Definition: ModError.F90:483
subroutine deregisterfunction(global)
Definition: ModError.F90:469