c xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx c c Flux smoothness based detection c c xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx subroutine cles_dcflag(ux,vx,dcflag,ncomps,nvars, $ ixlo,ixhi, nx, dx, direction,extra) implicit none integer ncomps, nvars integer ixlo, ixhi integer direction, extra(*) INTEGER nx integer dcflag(1:nx+1,1) double precision ux(ncomps,ixlo:ixhi) double precision vx(nvars,ixlo:ixhi) DOUBLE PRECISION dx integer i, m, slow, shigh, loop, size, enoOrder double precision smoothness, dmdx, maxv(4), minv(4), variance(4) integer span double precision TriggerValue, smooth_eps TriggerValue = 0.0d0 smooth_eps = 0.001 span = 3 enoOrder = extra(1) size = enoOrder-1 ! ---- test criteria at each point DO i = 1, nx+1, 1 DO loop =1,4 do m=i-size, i+size-1 dmdx = abs(fx(loop,m+1)-fx(loop,m-1)) if ( m .eq. i-size ) then maxv(loop) = dmdx minv(loop) = dmdx else maxv(loop) = MAX(maxv(loop), dmdx) minv(loop) = MIN(minv(loop), dmdx) endif enddo variance(loop) = maxv(loop)/(minv(loop)+smooth_eps) if ( loop .eq. 1 ) then smoothness = variance(loop) else smoothness = MAX(smoothness, variance(loop)) endif END DO IF ( smoothness .gt. TriggerValue ) THEN ! we are at a shock ! ---- we will use WENO here slow = max(i-span,1) shigh = min(i+span,nx+1) do m=slow,shigh dcflag(m,direction) = 1 enddo END IF END DO return end