Functions
	exact_b (ja, rr)

	run (n, ratio, ja, x0p=0.10, xep=0.25, alpha_l=1.0, cfl=0.2, fo=0.5, nt=4, kratio=None, rcpratio=None, consistent=True, plane=True, quad=False, muscl=False, verbose=True, energy_order=1, deposit=None)

	main ()

Detailed Description

WO-P23 rung P2 — the SUCKING INTERFACE (Welch & Wilson, JCP 160:662, 2000).

A planar interface with SATURATED VAPOUR on one side and SUPERHEATED LIQUID on the other. All the
latent heat comes from the liquid, so this is the first rung on which the LIQUID half of the
one-sided gradient fit, the per-phase k(C) / rho c_p(C) closures and the consistent rho c_p T
transport are all live at once (rung P1's liquid was saturated and pinned, so none of them entered).

GEOMETRY.  Vapour occupies x < X(t) against a wall held at T_sat; liquid occupies x > X(t) and is
superheated to T_inf far away.  The interface moves INTO the liquid and the liquid is pushed AWAY
from it at u_l = mdot (1/rho_v - 1/rho_l), so the domain needs an OUTLET at the liquid end: the net
vapour production is a real volume flux that has to leave.  (WO-P01's closed-column trick with a
balancing sink cannot be used here — the sink removes VOLUME from cells whose colour is 1, so the
colour drops there, and over a run this long the required removal is many cells' worth.)

WHY THE OUTLET IS SOUND AT HIGH DENSITY RATIO.  `applyBoundaryOpenness` re-imposes the literal 1.0
at a boundary face instead of the variable-density coefficient `open * rho0/rho_f` (WO-R2 item 1).
Here the outlet sits in PURE LIQUID, so choosing the reference density `set_rho(rho_l)` makes
`rho0/rho_f == 1` exactly at that face and the overwrite is the correct coefficient — the defect is
side-stepped by construction rather than tolerated.  Reported with every run.

THE EXACT SOLUTION (derived in the docstring of `exact_b`, checked against Welch & Wilson):

    X(t) = 2 b sqrt(alpha_l t) ,     b = Ja exp(-g^2) / (sqrt(pi) erfc(g)) ,  g = b rho_v/rho_l
    Ja   = rho_l c_pl (T_inf - T_sat) / (rho_v h_lv)
    T(x,t) = T_inf - (T_inf - T_sat) erfc(s)/erfc(g),  s = x/(2 sqrt(alpha_l t)) - b(1 - rho_v/rho_l)

Usage:
    PYTHONPATH=$PWD/build_cuda python tests/study/vof_sucking.py [--ratio 10] [--ns 64,128,256]

Function Documentation

◆ exact_b()

vof_sucking.exact_b	(	ja,
		rr
	)

Solve b = Ja exp(-g^2)/(sqrt(pi) erfc(g)), g = b*rr,  rr = rho_v/rho_l.

Derivation.  The vapour is isothermal at T_sat (it is bounded by a wall at T_sat), so it is at
rest and mdot = rho_v Xdot.  The liquid moves at the uniform u_l = Xdot (1 - rho_v/rho_l), so in
the Galilean frame x' = x - integral(u_l) the liquid temperature obeys the plain heat equation
and is an erfc in x'/(2 sqrt(alpha_l t)).  Matching T(X) = T_sat fixes the amplitude and the
Stefan condition mdot h_lv = k_l dT/dx|_X gives the fixed point above.

Definition at line 40 of file vof_sucking.py.

Referenced by run().

◆ run()

vof_sucking.run	(	n,
		ratio,
		ja,
		x0p = `0.10`,
		xep = `0.25`,
		alpha_l = `1.0`,
		cfl = `0.2`,
		fo = `0.5`,
		nt = `4`,
		kratio = `None`,
		rcpratio = `None`,
		consistent = `True`,
		plane = `True`,
		quad = `False`,
		muscl = `False`,
		verbose = `True`,
		energy_order = `1`,
		deposit = `None`
	)

One resolution.  The PHYSICAL problem is fixed (domain length 1, alpha_l = 1, the vapour
layer growing from x0p to xep of the domain) and only h = 1/n changes.  `nt` is the transverse
extent: the problem is 1-D so 4 is enough physically, but it decides how far the pressure
multigrid can coarsen (see the findings).

Definition at line 56 of file vof_sucking.py.

References exact_b().

Referenced by main().

◆ main()

vof_sucking.main ( )

Definition at line 178 of file vof_sucking.py.