47 INTEGER,
DIMENSION(:),
ALLOCATABLE :: cbtype
48 REAL,
DIMENSION(:,:,:),
ALLOCATABLE :: rscale, & !< radial scaling constants
50 REAL,
DIMENSION(:,:,:),
POINTER :: radius
92 INTEGER,
INTENT(IN) :: dir
100 ALLOCATE(this%cbtype(physics%VNUM),stat=err)
102 CALL this%Error(
"InitBoundary_custom",
"Unable to allocate memory.")
120 REAL,
INTENT(IN) :: time
125 SELECT CASE(this%GetDirection())
131 SELECT CASE(this%cbtype(m))
133 DO k=mesh%KMIN,mesh%KMAX
134 DO j=mesh%JMIN,mesh%JMAX
137 pvar%data4d(mesh%IMIN-i,j,k,m) = pvar%data4d(mesh%IMIN,j,k,m)
142 DO k=mesh%KMIN,mesh%KMAX
143 DO j=mesh%JMIN,mesh%JMAX
146 pvar%data4d(mesh%IMIN-i,j,k,m) = pvar%data4d(mesh%IMAX-i+1,j,k,m)
151 DO k=mesh%KMIN,mesh%KMAX
152 DO j=mesh%JMIN,mesh%JMAX
155 pvar%data4d(mesh%IMIN-i,j,k,m) = pvar%data4d(mesh%IMIN+i-1,j,k,m)
160 DO k=mesh%KMIN,mesh%KMAX
161 DO j=mesh%JMIN,mesh%JMAX
164 pvar%data4d(mesh%IMIN-i,j,k,m) = -pvar%data4d(mesh%IMIN+i-1,j,k,m)
169 DO k=mesh%KMIN,mesh%KMAX
170 DO j=mesh%JMIN,mesh%JMAX
173 pvar%data4d(mesh%IMIN-i,j,k,m) = (i+1)*pvar%data4d(mesh%IMIN,j,k,m) &
174 - i*pvar%data4d(mesh%IMIN+1,j,k,m)
179 DO k=mesh%KMIN,mesh%KMAX
180 DO j=mesh%JMIN,mesh%JMAX
183 pvar%data4d(mesh%IMIN-i,j,k,m) = this%data(i,j,k,m)
188 DO k=mesh%KMIN,mesh%KMAX
189 DO j=mesh%JMIN,mesh%JMAX
192 pvar%data4d(mesh%IMIN-i,j,k,m) = pvar%data4d(mesh%IMIN,j,k,m) &
193 * abs(pvar%data4d(mesh%IMIN,j,k,m) / pvar%data4d(mesh%IMIN+1,j,k,m))**i
198 DO k=mesh%KMIN,mesh%KMAX
199 DO j=mesh%JMIN,mesh%JMAX
200 IF (pvar%data4d(mesh%IMIN,j,k,m).GE.0.0 )
THEN
204 pvar%data4d(mesh%IMIN-i,j,k,m) = -pvar%data4d(mesh%IMIN+i-1,j,k,m)
210 pvar%data4d(mesh%IMIN-i,j,k,m) = pvar%data4d(mesh%IMIN,j,k,m)
216 DO k=mesh%KMIN,mesh%KMAX
217 DO j=mesh%JMIN,mesh%JMAX
220 pvar%data4d(mesh%IMIN-i,j,k,m) = (pvar%data4d(mesh%IMIN,j,k,m) &
221 + mesh%radius%bcenter(mesh%IMIN,j,k)*mesh%OMEGA) * this%invRscale(i,j,k) &
222 - mesh%radius%bcenter(mesh%IMIN-i,j,k)*mesh%OMEGA
227 DO k=mesh%KMIN,mesh%KMAX
228 DO j=mesh%JMIN,mesh%JMAX
231 pvar%data4d(mesh%IMIN-i,j,k,m) = (pvar%data4d(mesh%IMIN,j,k,m) &
232 + mesh%radius%bcenter(mesh%IMIN,j,k)**2*mesh%OMEGA) * this%Rscale(i,j,k) &
233 - mesh%radius%bcenter(mesh%IMIN-i,j,k)**2*mesh%OMEGA
241 this%err = ior(
custom,int(z
'0100'))
249 SELECT CASE(this%cbtype(m))
251 DO k=mesh%KMIN,mesh%KMAX
252 DO j=mesh%JMIN,mesh%JMAX
255 pvar%data4d(mesh%IMAX+i,j,k,m) = pvar%data4d(mesh%IMAX,j,k,m)
260 DO k=mesh%KMIN,mesh%KMAX
261 DO j=mesh%JMIN,mesh%JMAX
264 pvar%data4d(mesh%IMAX+i,j,k,m) = pvar%data4d(mesh%IMIN+i-1,j,k,m)
269 DO k=mesh%KMIN,mesh%KMAX
270 DO j=mesh%JMIN,mesh%JMAX
273 pvar%data4d(mesh%IMAX+i,j,k,m) = pvar%data4d(mesh%IMAX-i+1,j,k,m)
278 DO k=mesh%KMIN,mesh%KMAX
279 DO j=mesh%JMIN,mesh%JMAX
282 pvar%data4d(mesh%IMAX+i,j,k,m) = -pvar%data4d(mesh%IMAX-i+1,j,k,m)
287 DO k=mesh%KMIN,mesh%KMAX
288 DO j=mesh%JMIN,mesh%JMAX
291 pvar%data4d(mesh%IMAX+i,j,k,m) = (i+1)*pvar%data4d(mesh%IMAX,j,k,m) &
292 - i*pvar%data4d(mesh%IMAX-1,j,k,m)
297 DO k=mesh%KMIN,mesh%KMAX
298 DO j=mesh%JMIN,mesh%JMAX
301 pvar%data4d(mesh%IMAX+i,j,k,m) = this%data(i,j,k,m)
306 DO k=mesh%KMIN,mesh%KMAX
307 DO j=mesh%JMIN,mesh%JMAX
310 pvar%data4d(mesh%IMAX+i,j,k,m) = pvar%data4d(mesh%IMAX,j,k,m) &
311 * abs(pvar%data4d(mesh%IMAX,j,k,m) / pvar%data4d(mesh%IMAX-1,j,k,m))**i
316 DO k=mesh%KMIN,mesh%KMAX
317 DO j=mesh%JMIN,mesh%JMAX
319 IF (pvar%data4d(mesh%IMAX,j,k,m).LE.0.0)
THEN
321 pvar%data4d(mesh%IMAX+i,j,k,m) = -pvar%data4d(mesh%IMAX-i+1,j,k,m)
325 pvar%data4d(mesh%IMAX+i,j,k,m) = pvar%data4d(mesh%IMAX,j,k,m)
331 DO k=mesh%KMIN,mesh%KMAX
332 DO j=mesh%JMIN,mesh%JMAX
335 pvar%data4d(mesh%IMAX+i,j,k,m) = (pvar%data4d(mesh%IMAX,j,k,m) &
336 + mesh%radius%bcenter(mesh%IMAX,j,k)*mesh%OMEGA) * this%invRscale(i,j,k) &
337 - mesh%radius%bcenter(mesh%IMAX+i,j,k)*mesh%OMEGA
342 DO k=mesh%KMIN,mesh%KMAX
343 DO j=mesh%JMIN,mesh%JMAX
346 pvar%data4d(mesh%IMAX+i,j,k,m) = (pvar%data4d(mesh%IMAX,j,k,m) &
347 + mesh%radius%bcenter(mesh%IMAX,j,k)**2*mesh%OMEGA) * this%Rscale(i,j,k) &
348 - mesh%radius%bcenter(mesh%IMAX+i,j,k)**2*mesh%OMEGA
353 this%err = ior(
custom,int(z
'0200'))
358 SELECT CASE(this%cbtype(m))
362 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN-j,mesh%KMIN:mesh%KMAX,m) &
363 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN,mesh%KMIN:mesh%KMAX,m)
368 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN-j,mesh%KMIN:mesh%KMAX,m) &
369 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX-j+1,mesh%KMIN:mesh%KMAX,m)
374 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN-j,mesh%KMIN:mesh%KMAX,m) &
375 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN+j-1,mesh%KMIN:mesh%KMAX,m)
380 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN-j,mesh%KMIN:mesh%KMAX,m) &
381 = -pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN+j-1,mesh%KMIN:mesh%KMAX,m)
386 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN-j,mesh%KMIN:mesh%KMAX,m) &
387 = (j+1)*pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN,mesh%KMIN:mesh%KMAX,m) &
388 - j*pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN+1,mesh%KMIN:mesh%KMAX,m)
393 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN-j,mesh%KMIN:mesh%KMAX,m) &
394 = this%data(mesh%IMIN:mesh%IMAX,j,mesh%KMIN:mesh%KMAX,m)
399 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN-j,mesh%KMIN:mesh%KMAX,m) &
400 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN,mesh%KMIN:mesh%KMAX,m) &
401 * abs(pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN,mesh%KMIN:mesh%KMAX,m) &
402 / pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN+1,mesh%KMIN:mesh%KMAX,m))**j
416 WHERE (pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN,mesh%KMIN:mesh%KMAX,m).GE.0.0 )
417 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JGMIN,mesh%KMIN:mesh%KMAX,m) &
418 = -pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN+1,mesh%KMIN:mesh%KMAX,m)
419 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JGMIN+1,mesh%KMIN:mesh%KMAX,m) &
420 = -pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN,mesh%KMIN:mesh%KMAX,m)
422 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JGMIN,mesh%KMIN:mesh%KMAX,m) &
423 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN,mesh%KMIN:mesh%KMAX,m)
424 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JGMIN+1,mesh%KMIN:mesh%KMAX,m) &
425 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN,mesh%KMIN:mesh%KMAX,m)
430 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN-j,mesh%KMIN:mesh%KMAX,m) &
431 = (pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN,mesh%KMIN:mesh%KMAX,m) &
432 + mesh%radius%bcenter(mesh%IMIN:mesh%IMAX,mesh%JMIN,mesh%KMIN:mesh%KMAX)*mesh%OMEGA)&
433 * this%invRscale(mesh%IMIN:mesh%IMAX,j,mesh%KMIN:mesh%KMAX) &
434 - mesh%radius%bcenter(mesh%IMIN:mesh%IMAX,mesh%JMIN-j,mesh%KMIN:mesh%KMAX)*mesh%OMEGA
439 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN-j,mesh%KMIN:mesh%KMAX,m) &
440 = (pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN,mesh%KMIN:mesh%KMAX,m) &
441 + mesh%radius%bcenter(mesh%IMIN:mesh%IMAX,mesh%JMIN,mesh%KMIN:mesh%KMAX)**2*mesh%OMEGA)&
442 * this%Rscale(mesh%IMIN:mesh%IMAX,j,mesh%KMIN:mesh%KMAX) &
443 - mesh%radius%bcenter(mesh%IMIN:mesh%IMAX,mesh%JMIN-j,mesh%KMIN:mesh%KMAX)**2*mesh%OMEGA
446 this%err = ior(
custom,int(z
'0300'))
451 SELECT CASE(this%cbtype(m))
455 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX+j,mesh%KMIN:mesh%KMAX,m) &
456 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX,mesh%KMIN:mesh%KMAX,m)
461 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX+j,mesh%KMIN:mesh%KMAX,m) &
462 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN+j-1,mesh%KMIN:mesh%KMAX,m)
467 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX+j,mesh%KMIN:mesh%KMAX,m) &
468 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX-j+1,mesh%KMIN:mesh%KMAX,m)
473 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX+j,mesh%KMIN:mesh%KMAX,m) &
474 = -pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX-j+1,mesh%KMIN:mesh%KMAX,m)
479 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX+j,mesh%KMIN:mesh%KMAX,m) &
480 = (j+1)*pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX,mesh%KMIN:mesh%KMAX,m) &
481 - j*pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX-1,mesh%KMIN:mesh%KMAX,m)
486 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX+j,mesh%KMIN:mesh%KMAX,m) &
487 = this%data(mesh%IMIN:mesh%IMAX,j,mesh%KMIN:mesh%KMAX,m)
492 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX+j,mesh%KMIN:mesh%KMAX,m) &
493 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX,mesh%KMIN:mesh%KMAX,m) &
494 * abs(pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX,mesh%KMIN:mesh%KMAX,m) &
495 / pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX-1,mesh%KMIN:mesh%KMAX,m))**j
509 WHERE (pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX,mesh%KMIN:mesh%KMAX,m).LE.0.0 )
510 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JGMAX,mesh%KMIN:mesh%KMAX,m) &
511 = -pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX-1,mesh%KMIN:mesh%KMAX,m)
512 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JGMAX-1,mesh%KMIN:mesh%KMAX,m) &
513 = -pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX,mesh%KMIN:mesh%KMAX,m)
515 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JGMAX,mesh%KMIN:mesh%KMAX,m) &
516 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX,mesh%KMIN:mesh%KMAX,m)
517 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JGMAX-1,mesh%KMIN:mesh%KMAX,m) &
518 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX,mesh%KMIN:mesh%KMAX,m)
523 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX+j,mesh%KMIN:mesh%KMAX,m) &
524 = (pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX,mesh%KMIN:mesh%KMAX,m) &
525 + mesh%radius%bcenter(mesh%IMIN:mesh%IMAX,mesh%JMAX,mesh%KMIN:mesh%KMAX)*mesh%OMEGA)&
526 * this%invRscale(mesh%IMIN:mesh%IMAX,j,mesh%KMIN:mesh%KMAX) &
527 - mesh%radius%bcenter(mesh%IMIN:mesh%IMAX,mesh%JMAX+j,mesh%KMIN:mesh%KMAX)*mesh%OMEGA
532 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX+j,mesh%KMIN:mesh%KMAX,m) &
533 = (pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMAX,mesh%KMIN:mesh%KMAX,m) &
534 + mesh%radius%bcenter(mesh%IMIN:mesh%IMAX,mesh%JMAX,mesh%KMIN:mesh%KMAX)**2*mesh%OMEGA)&
535 * this%Rscale(mesh%IMIN:mesh%IMAX,j,mesh%KMIN:mesh%KMAX) &
536 - mesh%radius%bcenter(mesh%IMIN:mesh%IMAX,mesh%JMAX+j,mesh%KMIN:mesh%KMAX)**2*mesh%OMEGA
539 this%err = ior(
custom,int(z
'0400'))
544 SELECT CASE(this%cbtype(m))
548 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN-k,m) &
549 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN,m)
554 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN-k,m) &
555 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX-k+1,m)
560 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN-k,m) &
561 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN+k-1,m)
566 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN-k,m) &
567 = -pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN+k-1,m)
572 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN-k,m) &
573 = (k+1)*pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN,m) &
574 - k*pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN+1,m)
579 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN-k,m) &
580 = this%data(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,k,m)
585 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN-k,m) &
586 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN,m) &
587 * abs(pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN,m) &
588 / pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN+1,m))**k
604 WHERE (pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN,m).GE.0.0 )
605 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KGMIN,m) &
606 = -pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN+1,m)
607 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KGMIN+1,m) &
608 = -pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN,m)
610 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KGMIN,m) &
611 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN,m)
612 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KGMIN+1,m) &
613 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN,m)
619 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN-k,m) &
620 = (pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN,m) &
621 + mesh%radius%bcenter(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN)*mesh%OMEGA)&
622 * this%invRscale(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,k) &
623 - mesh%radius%bcenter(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN-k)*mesh%OMEGA
628 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN-k,m) &
629 = (pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN,m) &
630 + mesh%radius%bcenter(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN)**2*mesh%OMEGA)&
631 * this%invRscale(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,k) &
632 - mesh%radius%bcenter(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN-k)**2*mesh%OMEGA
635 this%err = ior(
custom,int(z
'0500'))
640 SELECT CASE(this%cbtype(m))
644 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX+k,m) &
645 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX,m)
650 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX+k,m) &
651 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN+k-1,m)
656 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX+k,m) &
657 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX-k+1,m)
662 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX+k,m) &
663 = -pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX-k+1,m)
668 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX+k,m) &
669 = (k+1)*pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX,m) &
670 - k*pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX-1,m)
675 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX+k,m) &
676 = this%data(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,k,m)
681 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX+k,m) &
682 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX,m) &
683 * abs(pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX,m) &
684 / pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX-1,m))**k
699 WHERE (pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX,m).LE.0.0)
700 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KGMAX,m) &
701 = -pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX-1,m)
702 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KGMAX-1,m) &
703 = -pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX,m)
705 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KGMAX,m) &
706 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX,m)
707 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KGMAX-1,m) &
708 = pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX,m)
714 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX+k,m) &
715 = (pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX,m) &
716 + mesh%radius%bcenter(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX)*mesh%OMEGA)&
717 * this%invRscale(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,k) &
718 - mesh%radius%bcenter(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX+k)*mesh%OMEGA
723 pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX+k,m) &
724 = (pvar%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX,m) &
725 + mesh%radius%bcenter(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX)**2*mesh%OMEGA)&
726 * this%invRscale(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,k) &
727 - mesh%radius%bcenter(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMAX+k)**2*mesh%OMEGA
730 this%err = ior(
custom,int(z
'0600'))
742 INTEGER,
DIMENSION(1:Physics%VNUM) :: cbtype
743 REAL,
DIMENSION(:,:,:),
POINTER,
OPTIONAL :: kepler_radius
745 INTEGER :: i,j,k,err = 0
747 this%cbtype(:) = cbtype(:)
749 SELECT CASE(this%GetDirection())
751 ALLOCATE(this%data(mesh%GINUM,mesh%JMIN:mesh%JMAX,mesh%KMIN:mesh%KMAX,physics%VNUM), &
754 ALLOCATE(this%data(mesh%IMIN:mesh%IMAX,mesh%GJNUM,mesh%KMIN:mesh%KMAX,physics%VNUM), &
757 ALLOCATE(this%data(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%GKNUM,physics%VNUM), &
760 IF (err.NE.0)
CALL this%Error(
"SetCustomBoundaries",
"Unable to allocate memory.")
761 this%data(:,:,:,:) = 0.0
765 SELECT CASE(this%GetDirection())
767 IF (mesh%GINUM.NE.2)
CALL this%Error(
"SetCustomBoundaries",
"GINUM must be 2 for outflow boundaries")
769 IF (mesh%GJNUM.NE.2)
CALL this%Error(
"SetCustomBoundaries",
"GJNUM must be 2 for outflow boundaries")
771 IF (mesh%GKNUM.NE.2)
CALL this%Error(
"SetCustomBoundaries",
"GKNUM must be 2 for outflow boundaries")
776 SELECT CASE(this%GetDirection())
778 ALLOCATE(this%invRscale(mesh%GINUM,mesh%JMIN:mesh%JMAX,mesh%KMIN:mesh%KMAX), &
781 ALLOCATE(this%invRscale(mesh%IMIN:mesh%IMAX,mesh%GJNUM,mesh%KMIN:mesh%KMAX), &
784 ALLOCATE(this%invRscale(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%GKNUM), &
787 IF (err.NE.0)
CALL this%Error(
"SetCustomBoundaries",
"Unable to allocate memory.")
789 IF (
PRESENT(kepler_radius))
THEN
790 this%radius => kepler_radius
791 IF (mesh%OMEGA.GT.tiny(mesh%OMEGA)) &
792 CALL this%Warning(
"SetCustomBoundaries", &
793 "user supplied radius and rotating frame may yield unexpected/wrong results")
795 this%radius => mesh%radius%bcenter
798 SELECT CASE(this%GetDirection())
800 FORALL (i=1:mesh%GINUM,j=mesh%JMIN:mesh%JMAX,k=mesh%KMIN:mesh%KMAX) &
801 this%invRscale(i,j,k) = this%radius(mesh%IMIN,j,k) / this%radius(mesh%IMIN-i,j,k)
803 FORALL (i=1:mesh%GINUM,j=mesh%JMIN:mesh%JMAX,k=mesh%KMIN:mesh%KMAX) &
804 this%invRscale(i,j,k) = this%radius(mesh%IMAX,j,k) / this%radius(mesh%IMAX+i,j,k)
806 FORALL (i=mesh%IMIN:mesh%IMAX,j=1:mesh%GjNUM,k=mesh%KMIN:mesh%KMAX) &
807 this%invRscale(i,j,k) = this%radius(i,mesh%JMIN,k) / this%radius(i,mesh%JMIN-j,k)
809 FORALL (i=mesh%IMIN:mesh%IMAX,j=1:mesh%GjNUM,k=mesh%KMIN:mesh%KMAX) &
810 this%invRscale(i,j,k) = this%radius(i,mesh%JMAX,k) / this%radius(i,mesh%JMAX+j,k)
812 FORALL (i=mesh%IMIN:mesh%IMAX,j=mesh%JMIN:mesh%JMAX,k=1:mesh%GKNUM) &
813 this%invRscale(i,j,k) = this%radius(i,j,mesh%KMIN) / this%radius(i,j,mesh%KMIN-k)
815 FORALL (i=mesh%IMIN:mesh%IMAX,j=mesh%JMIN:mesh%JMAX,k=1:mesh%GKNUM) &
816 this%invRscale(i,j,k) = this%radius(i,j,mesh%KMAX) / this%radius(i,j,mesh%KMAX+k)
818 this%invRscale(:,:,:) = sqrt(this%invRscale(:,:,:))
822 SELECT CASE(this%GetDirection())
824 ALLOCATE(this%Rscale(mesh%GINUM,mesh%JMIN:mesh%JMAX,mesh%KMIN:mesh%KMAX), &
827 ALLOCATE(this%Rscale(mesh%IMIN:mesh%IMAX,mesh%GJNUM,mesh%KMIN:mesh%KMAX), &
830 ALLOCATE(this%Rscale(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%GKNUM), &
834 IF (
PRESENT(kepler_radius))
THEN
835 this%radius => kepler_radius
836 IF (mesh%OMEGA.GT.tiny(mesh%OMEGA)) &
837 CALL this%Warning(
"SetCustomBoundaries", &
838 "user supplied radius and rotating frame may yield unexpected/wrong results")
840 this%radius => mesh%radius%bcenter
843 IF (err.NE.0)
CALL this%Error(
"SetCustomBoundaries",
"Unable to allocate memory.")
844 SELECT CASE(this%GetDirection())
846 FORALL (i=1:mesh%GINUM,j=mesh%JMIN:mesh%JMAX,k=mesh%KMIN:mesh%KMAX) &
847 this%Rscale(i,j,k) = this%radius(mesh%IMIN-i,j,k) / this%radius(mesh%IMIN,j,k)
849 FORALL (i=1:mesh%GINUM,j=mesh%JMIN:mesh%JMAX,k=mesh%KMIN:mesh%KMAX) &
850 this%Rscale(i,j,k) = this%radius(mesh%IMAX+i,j,k) / this%radius(mesh%IMAX,j,k)
852 FORALL (i=mesh%IMIN:mesh%IMAX,j=1:mesh%GjNUM,k=mesh%KMIN:mesh%KMAX) &
853 this%Rscale(i,j,k) = this%radius(i,mesh%JMIN-j,k) / this%radius(i,mesh%JMIN,k)
855 FORALL (i=mesh%IMIN:mesh%IMAX,j=1:mesh%GjNUM,k=mesh%KMIN:mesh%KMAX) &
856 this%Rscale(i,j,k) = this%radius(i,mesh%JMAX+j,k) / this%radius(i,mesh%JMAX,k)
858 FORALL (i=mesh%IMIN:mesh%IMAX,j=mesh%JMIN:mesh%JMAX,k=1:mesh%GKNUM) &
859 this%Rscale(i,j,k) = this%radius(i,j,mesh%KMIN-k) / this%radius(i,j,mesh%KMIN)
861 FORALL (i=mesh%IMIN:mesh%IMAX,j=mesh%JMIN:mesh%JMAX,k=1:mesh%GKNUM) &
862 this%Rscale(i,j,k) = this%radius(i,j,mesh%KMAX+k) / this%radius(i,j,mesh%KMAX)
864 this%Rscale(:,:,:) = sqrt(this%Rscale(:,:,:))
873 IF (
ALLOCATED(this%data))
DEALLOCATE(this%data)
874 IF (
ALLOCATED(this%cbtype))
DEALLOCATE(this%cbtype)
875 IF (
ALLOCATED(this%Rscale))
DEALLOCATE(this%Rscale)
876 IF (
ALLOCATED(this%invRscale))
DEALLOCATE(this%invRscale)
878 CALL this%Finalize_base()
integer, parameter custom
user defined
Boundary module for custom conditions.
@, public custom_extrapol
subroutine setboundarydata(this, Mesh, Physics, time, pvar)
Applies the custom boundary conditions.
subroutine finalize(this)
character(len=32), parameter boundcond_name
subroutine setcustomboundaries(this, Mesh, Physics, cbtype, kepler_radius)
@, public custom_angkepler
subroutine initboundary_custom(this, Mesh, Physics, dir, config)
Constructor for custom boundary conditions.
@, public custom_logexpol
Boundary module for fixed in/outflow conditions.
Dictionary for generic data types.
derived class for compound of mesh arrays
integer, parameter east
named constant for eastern boundary
integer, parameter bottom
named constant for bottom boundary
integer, parameter south
named constant for southern boundary
integer, parameter top
named constant for top boundary
integer, parameter north
named constant for northern boundary
integer, parameter west
named constant for western boundary