c c c ===================================================== subroutine bc3(maxmx,maxmy,maxmz,meqn,mbc,mx,my,mz,xlower, & ylower,zlower,dx,dy,dz,q,maux,aux,t,dt,mthbc) c ===================================================== c c # Standard boundary condition choices for claw3 c c # At each boundary k = 1 (xlower), 2 (xupper), c # 3 (ylower), 4 (yupper), c # 5 (zlower), 6 (zupper): c # mthbc(k) = 0 for user-supplied BC's (must be inserted!) c # = 1 for zero-order extrapolation c # = 2 for periodic boundary coniditions c # = 3 for solid walls, assuming this can be implemented c # by reflecting the data about the boundary and then c # negating the 2'nd (for k=1,2) or 3'rd (for k=3,4) c # or 4'th (for k=5,6) component of q. c ------------------------------------------------ c c # Extend the data from the interior cells (1:mx, 1:my, 1:mz) c # to a layer of mbc ghost cells outside the region. c implicit double precision (a-h,o-z) dimension q(1-mbc:maxmx+mbc, 1-mbc:maxmy+mbc, & 1-mbc:maxmz+mbc, meqn) dimension aux(1-mbc:maxmx+mbc, 1-mbc:maxmy+mbc, & 1-mbc:maxmz+mbc, maux) dimension mthbc(6) c c c------------------------------------------------------- c # left boundary (xlower): c------------------------------------------------------- go to (100,110,120,130) mthbc(1)+1 c 100 continue c # user-specified boundary conditions go here in place of error output write(6,*) '*** ERROR *** mthbc(1)=0 and no BCs specified in bc3' stop go to 199 c 110 continue c # zero-order extrapolation: do 115 m=1,meqn do 115 ibc=1,mbc do 115 j = 1-mbc, my+mbc do 115 k = 1-mbc, mz+mbc q(1-ibc,j,k,m) = q(1,j,k,m) 115 continue go to 199 120 continue c # periodic: do 125 m=1,meqn do 125 ibc=1,mbc do 125 j = 1-mbc, my+mbc do 125 k = 1-mbc, mz+mbc q(1-ibc,j,k,m) = q(mx+1-ibc,j,k,m) 125 continue go to 199 130 continue c # solid wall (assumes 2'nd component is velocity or momentum in x): do 135 m=1,meqn do 135 ibc=1,mbc do 135 j = 1-mbc, my+mbc do 135 k = 1-mbc, mz+mbc q(1-ibc,j,k,m) = q(ibc,j,k,m) 135 continue c # negate the normal velocity: do 136 ibc=1,mbc do 136 j = 1-mbc, my+mbc do 136 k = 1-mbc, mz+mbc q(1-ibc,j,k,2) = -q(ibc,j,k,2) 136 continue go to 199 199 continue c c------------------------------------------------------- c # right boundary (xupper): c------------------------------------------------------- go to (200,210,220,230) mthbc(2)+1 c 200 continue c # user-specified boundary conditions go here in place of error output write(6,*) '*** ERROR *** mthbc(2)=0 and no BCs specified in bc3' stop go to 299 210 continue c # zero-order extrapolation: do 215 m=1,meqn do 215 ibc=1,mbc do 215 j = 1-mbc, my+mbc do 215 k = 1-mbc, mz+mbc q(mx+ibc,j,k,m) = q(mx,j,k,m) 215 continue go to 299 220 continue c # periodic: do 225 m=1,meqn do 225 ibc=1,mbc do 225 j = 1-mbc, my+mbc do 225 k = 1-mbc, mz+mbc q(mx+ibc,j,k,m) = q(ibc,j,k,m) 225 continue go to 299 230 continue c # solid wall (assumes 2'nd component is velocity or momentum in x): do 235 m=1,meqn do 235 ibc=1,mbc do 235 j = 1-mbc, my+mbc do 235 k = 1-mbc, mz+mbc q(mx+ibc,j,k,m) = q(mx+1-ibc,j,k,m) 235 continue c # negate the normal velocity: do 236 ibc=1,mbc do 236 j = 1-mbc, my+mbc do 236 k = 1-mbc, mz+mbc q(mx+ibc,j,k,2) = -q(mx+1-ibc,j,k,2) 236 continue go to 299 299 continue c c------------------------------------------------------- c # bottom boundary (ylower): c------------------------------------------------------- go to (300,310,320,330) mthbc(3)+1 c 300 continue c # user-specified boundary conditions go here in place of error output write(6,*) '*** ERROR *** mthbc(3)=0 and no BCs specified in bc3' stop go to 399 c 310 continue c # zero-order extrapolation: do 315 m=1,meqn do 315 jbc=1,mbc do 315 i = 1-mbc, mx+mbc do 315 k = 1-mbc, mz+mbc q(i,1-jbc,k,m) = q(i,1,k,m) 315 continue go to 399 320 continue c # periodic: do 325 m=1,meqn do 325 jbc=1,mbc do 325 i = 1-mbc, mx+mbc do 325 k = 1-mbc, mz+mbc q(i,1-jbc,k,m) = q(i,my+1-jbc,k,m) 325 continue go to 399 330 continue c # solid wall (assumes 3'rd component is velocity or momentum in y): do 335 m=1,meqn do 335 jbc=1,mbc do 335 i = 1-mbc, mx+mbc do 335 k = 1-mbc, mz+mbc q(i,1-jbc,k,m) = q(i,jbc,k,m) 335 continue c # negate the normal velocity: do 336 jbc=1,mbc do 336 i = 1-mbc, mx+mbc do 336 k = 1-mbc, mz+mbc q(i,1-jbc,k,3) = -q(i,jbc,k,3) 336 continue go to 399 399 continue c c------------------------------------------------------- c # top boundary (yupper): c------------------------------------------------------- go to (400,410,420,430) mthbc(4)+1 c 400 continue c # user-specified boundary conditions go here in place of error output write(6,*) '*** ERROR *** mthbc(4)=0 and no BCs specified in bc3' stop go to 499 410 continue c # zero-order extrapolation: do 415 m=1,meqn do 415 jbc=1,mbc do 415 i = 1-mbc, mx+mbc do 415 k = 1-mbc, mz+mbc q(i,my+jbc,k,m) = q(i,my,k,m) 415 continue go to 499 420 continue c # periodic: do 425 m=1,meqn do 425 jbc=1,mbc do 425 i = 1-mbc, mx+mbc do 425 k = 1-mbc, mz+mbc q(i,my+jbc,k,m) = q(i,jbc,k,m) 425 continue go to 499 430 continue c # solid wall (assumes 3'rd component is velocity or momentum in y): do 435 m=1,meqn do 435 jbc=1,mbc do 435 i = 1-mbc, mx+mbc do 435 k = 1-mbc, mz+mbc q(i,my+jbc,k,m) = q(i,my+1-jbc,k,m) 435 continue c # negate the normal velocity: do 436 jbc=1,mbc do 436 i = 1-mbc, mx+mbc do 436 k = 1-mbc, mz+mbc q(i,my+jbc,k,3) = -q(i,my+1-jbc,k,3) 436 continue go to 499 499 continue c c------------------------------------------------------- c # boundary (zlower): c------------------------------------------------------- go to (500,510,520,530) mthbc(5)+1 c 500 continue c # user-specified boundary conditions go here in place of error output write(6,*) '*** ERROR *** mthbc(5)=0 and no BCs specified in bc3' stop go to 599 c 510 continue c # zero-order extrapolation: do 515 m=1,meqn do 515 kbc=1,mbc do 515 i = 1-mbc, mx+mbc do 515 j = 1-mbc, my+mbc q(i,j,1-kbc,m) = q(i,j,1,m) 515 continue go to 599 520 continue c # periodic: do 525 m=1,meqn do 525 kbc=1,mbc do 525 i = 1-mbc, mx+mbc do 525 j = 1-mbc, my+mbc q(i,j,1-kbc,m) = q(i,j,mz+1-kbc,m) 525 continue go to 599 530 continue c # solid wall (assumes 4'rd component is velocity or momentum in y): do 535 m=1,meqn do 535 kbc=1,mbc do 535 i = 1-mbc, mx+mbc do 535 j = 1-mbc, my+mbc q(i,j,1-kbc,m) = q(i,j,kbc,m) 535 continue c # negate the normal velocity: do 536 kbc=1,mbc do 536 i = 1-mbc, mx+mbc do 536 j = 1-mbc, my+mbc q(i,j,1-kbc,4) = -q(i,j,kbc,4) 536 continue go to 599 599 continue c c------------------------------------------------------- c # boundary (zupper): c------------------------------------------------------- go to (600,610,620,630) mthbc(6)+1 c 600 continue c # user-specified boundary conditions go here in place of error output c write(6,*) '*** ERROR *** mthbc(6)=0 and no BCs specified in bc3' c stop c # solid wall with piston in center do 605 m=1,meqn do 605 kbc=1,mbc do 605 i = 1-mbc, mx+mbc do 605 j = 1-mbc, my+mbc q(i,j,mz+kbc,m) = q(i,j,mz+1-kbc,m) 605 continue c # negate the normal velocity: do 606 kbc=1,mbc do 606 i = 1-mbc, mx+mbc do 606 j = 1-mbc, my+mbc if (t .lt. 0.1d0) then xi = xlower + (i-0.5d0)*dx yj = ylower + (j-0.5d0)*dy r2 = (xi**2 + yj**2) s = -dexp(-500.d0*r2)*dexp(-500*t**2) c write(27,271) t,i,j,r2,s c 271 format(d16.6,2i5,2d16.6) if (dabs(s).lt.1d-6) s = 0.d0 else s = 0.d0 endif q(i,j,mz+kbc,4) = 2.d0*s - q(i,j,mz+1-kbc,4) 606 continue go to 699 610 continue c # zero-order extrapolation: do 615 m=1,meqn do 615 kbc=1,mbc do 615 i = 1-mbc, mx+mbc do 615 j = 1-mbc, my+mbc q(i,j,mz+kbc,m) = q(i,j,mz,m) 615 continue go to 699 620 continue c # periodic: do 625 m=1,meqn do 625 kbc=1,mbc do 625 i = 1-mbc, mx+mbc do 625 j = 1-mbc, my+mbc q(i,j,mz+kbc,m) = q(i,j,kbc,m) 625 continue go to 699 630 continue c # solid wall (assumes 3'rd component is velocity or momentum in y): do 635 m=1,meqn do 635 kbc=1,mbc do 635 i = 1-mbc, mx+mbc do 635 j = 1-mbc, my+mbc q(i,j,mz+kbc,m) = q(i,j,mz+1-kbc,m) 635 continue c # negate the normal velocity: do 636 kbc=1,mbc do 636 i = 1-mbc, mx+mbc do 636 j = 1-mbc, my+mbc q(i,j,mz+kbc,4) = -q(i,j,mz+1-kbc,4) 636 continue go to 699 699 continue return end