純粋なVolumetricが重すぎて敗北したのでポスプロで実装したバージョン（Twigl の geeker(MRT)で動作します）

Undersea_volumetric_for_twigl_geekerMRT.frag

// https://www.shadertoy.com/view/4sfGzS
float hash(vec3 p)  // replace this by something better
{
    p  = fract( p*0.3183099+.1 );
    p *= 17.0;
    return fract( p.x*p.y*p.z*(p.x+p.y+p.z) );
}

float noise( in vec3 x )
{
    vec3 i = floor(x);
    vec3 f = fract(x);
    f = f*f*(3.0-2.0*f);

    return mix(mix(mix( hash(i+vec3(0,0,0)),
                        hash(i+vec3(1,0,0)),f.x),
                   mix( hash(i+vec3(0,1,0)),
                        hash(i+vec3(1,1,0)),f.x),f.y),
               mix(mix( hash(i+vec3(0,0,1)),
                        hash(i+vec3(1,0,1)),f.x),
                   mix( hash(i+vec3(0,1,1)),
                        hash(i+vec3(1,1,1)),f.x),f.y),f.z);
}
float fbm( in vec3 x )
{
    const float scale  = 1.5;

    float a = 0.0;
    float b = 0.5;
    float f = 1.0;
    vec3  d = vec3(0.0);
    for( int i=0; i<8; i++ )
    {
        float n = noise(f*x*scale);
        a += b*n;             // accumulate values
        b *= 0.5;             // amplitude decrease
        f *= 1.8;             // frequency increase
    }

  return a;
}

// http://iquilezles.org/www/articles/smin/smin.htm
float smin( float a, float b, float k )
{
    float h = max(k-abs(a-b),0.0);
    return min(a, b) - h*h*0.25/k;
}
float smax( float a, float b, float k )
{
  float h = clamp( 0.5 + 0.5*(b-a)/k, 0.0, 1.0 );
  return mix( a, b, h ) + k*h*(1.0-h);
}

float map2( in vec3 p )
{
  float d = fbm(p * .1);
  d -= 0.3;
  d *= 9.;
  return d;
}

float map(vec3 p) {
  float d = map2(p);
  d = smin(d, p.y+12., 9.);
  d = smax(d, -length(p.xy)+8., 5.);
  p.y+=2.;
  p.x = abs(p.x);
  d = min(d, length(p.xy-vec2(2., -3)));
  return d;
}

vec3 light(vec3 lPos, vec3 p) {
    vec3 L = lPos-p;
    return (50000.0*vec3(0.3, 0.6, 1.2)) * (1.0/dot(L,L));
}

void pass1(out vec4 o) {
  vec2 p = (gl_FragCoord.xy*2. - r.xy) / max(r.x, r.y);
  p /= 1. - length(p) * 0.8;

  vec3 cPos = vec3(0,0,t*35.0);
  float dx = sin((cPos.z)*0.014)*3.0;
  float dy = sin((cPos.z)*0.02)*2.0+1.;
  cPos.xy -= 0.5*vec2(dx, dy);
  vec3 cDir = normalize(vec3(dx, dy, 40.0));
  vec3 cUp = vec3(0.0, 1.0, 0.0);
  vec3 cSide = cross(cUp, cDir);

  vec3 ray = normalize(cSide * p.x + cUp * p.y + cDir);

  vec3 lPos = vec3(0,50,50)+cPos;
  vec3 color = vec3(0);

  float depth = 0.;
  vec3 pos = vec3(0);
  for(int i=0; i<99; i++) {
    pos = cPos + ray * depth;
    float d = map(pos);
    color += 0.0007 * light(lPos, pos) ;
    depth += d;
    if (d<0.01 || depth>200.) break;
  }

  o = vec4(color, 1.);
}


// https://www.shadertoy.com/view/MdG3RD
void pass2(sampler2D b, out vec4 o) {
	// 2pass
  vec2 uv=gl_FragCoord.xy/r.xy;
  o = texture(b, uv);
  // return;
	float decay=0.84;
	float exposure=0.9;
	float density=0.8;
	float weight=0.5;

	vec2 tc = uv;
	vec2 lightPos = vec2(0.);
	vec2 deltaTexCoord = tc;

  const int NUM_SAMPLES = 55;
  deltaTexCoord = uv-normalize(vec2(1, 2));
	deltaTexCoord *= 1.0 / float(NUM_SAMPLES)  * density;

	float illuminationDecay = 1.0;
	vec4 color =texture(b, tc.xy)*0.305104;

	tc += deltaTexCoord * fract( sin(dot(uv.xy+fract(t),
																			 vec2(12.9898, 78.233)))* 43758.5453 );

	for(int i=0; i < NUM_SAMPLES; i++)
	{
			tc -= deltaTexCoord;
			vec4 sampleTex = texture(b, tc)*0.305104;
			sampleTex *= illuminationDecay * weight;
			color += sampleTex;
			illuminationDecay *= decay;
	}
	o = color*exposure;
}

void main() {
	pass1(o1);
	pass2(b1, o0);
}

Preview

https://twigl.app/

Tweet

https://twitter.com/kasari39/status/1261978601584394242