{"id":"3eee5183-6dea-4b4b-a59a-be19a86ee90a","arxiv_id":"2506.09077","paper_version":1,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":6.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":0,"one_line_summary":"The paper classifies all Riemann solutions for three-phase flow with heavy oil, where oil viscosity dominates, for left states on the water-gas edge and right states across nearly the whole saturation triangle.","lead":"This math paper works out the full list of possible wave patterns when water and gas are injected into a heavy-oil reservoir, for almost all possible starting states of the oil. It gives engineers a theoretical map of how oil, water, and gas fronts move, which helps in designing enhanced oil recovery and carbon storage projects.","discovery_kind":"extension","skeptic_critique":{"model":"deepseek-v4-flash","headline":"The classification relies on shock admissibility verified only for D(S)=identity; with a non-identity capillary pressure matrix the admissible Hugoniot segments, and hence the Riemann solution structure, can change.","rationale":"The paper's central contribution is a classification of Riemann solutions for left states on G-W and right states over most of the saturation triangle. Everything in that classification flows from which shocks on the Hugoniot loci are admissible. The admissibility criterion is the viscous profile criterion, but it is implemented with the simplifying choice D(S)=I (Section 3.1, Remark 3.1). This is not merely a technical convenience: the identity matrix makes the ODE (9) invariant along the secondary bifurcation segments G-D, W-E, O-B, and the paper relies on that invariance to decide nonlocal shock admissibility by relative position. For a general positive-definite capillary pressure matrix, that invariance fails, and the set of admissible shocks can change. Because the wave curves W_f(R) are built from admissible shock segments, the whole classification (Claims 4.1-4.30) is conditional on D=I. The abstract and conclusion do not state this condition, so the advertised scope ('three-phase flow in porous media') exceeds what is justified. The paper is transparent about the simplification, and the numerical simulations in Section 5 confirm the D=I model, so the work is valuable as a conditional classification; however, the central claim's generality is not established. A single computational experiment with a non-identity D on representative shocks would settle whether the classification is robust or model-specific. The reader's CONDITIONAL verdict is appropriate; no change is needed. I also note the abstract overstates the viscosity regime (only inequalities (5) are mentioned, while Section 2 adds a further restriction on the double contact locus), and L1_loc-stability is supported only numerically, but these are secondary to the D=I concern.","tokens_in":40262,"tokens_out":10353,"duration_ms":123584,"concrete_test":"Using the paper's representative viscosities (µ_w=1, µ_o=9.5, µ_g=0.45), pick R=(0.0402518,0.913397) in subregion Θ_1^a (Fig. 9(a)). For each shock on the segments [A1,R) and [A2,A3] declared admissible under D=I, numerically solve the traveling-wave ODE (9) with a non-identity positive-definite matrix, e.g., D=diag(1,0.45), and check whether a heteroclinic connection from M to R still exists. If any declared admissible segment loses its profile, or any previously inadmissible segment gains one, the wave curve W_f(R) changes and the classification of Section 4.1 depends on the identity-viscosity simplification.","verdict_should_be":"UNCHANGED","load_bearing_attack":"In Section 3.1 the authors state 'we assume D(S) to be the identity matrix and numerically verify the admissibility of shocks' (Eq. (9) with D=I); Remark 3.1 reiterates that every shock must satisfy the viscous profile criterion with identity viscosity. All admissible shock segments used in the wave curves, e.g. [A1,R) and [A2,A3] in region Θ_1 (Figs. 9-10), are selected by numerically checking existence of traveling-wave solutions of this D=I ODE. The identity choice is not innocuous: for non-identity positive-definite D, the ODE vector field changes and the segments G-D, W-E, O-B are no longer invariant, so the stated criterion that 'admissibility of nonlocal shocks is determined by the relative positions of M and N with respect to the corresponding segments' (Section 3.1) no longer applies. If a nonlocal shock that is admissible under D=I becomes inadmissible under a physical capillary pressure matrix, or vice versa, the backward fast wave curve W_f(R) changes, and with it the Riemann solution classifications in Claims 4.1-4.30. The abstract and conclusion present the results for three-phase flow in porous media without qualifying this D=I simplification, so the central claim as advertised depends on a modeling assumption that is not physically justified.","agreement_with_reader":"agree"},"referee_report":{"model":"deepseek-v4-flash","summary":"The paper studies the Riemann problem for a 2x2 system of conservation laws modeling three-phase immiscible flow in a porous medium with Corey-type quadratic relative permeabilities, in the heavy-oil regime where oil viscosity greatly exceeds water and gas viscosities. The left state L is taken on the G-W edge of the saturation triangle and the right state R covers nearly the whole saturation triangle except small regions near G-O and W-O. The authors use the wave curve method to construct solutions, dividing the right-state triangle into regions Θ, Ω, Γ (and subregions) and presenting the resulting wave-group structure in Claims 4.1-4.30. Shock admissibility is enforced through the viscous profile criterion, but with the viscosity matrix D(S) set to the identity and admissibility verified numerically. Three numerical simulations are compared with the analytical profiles. The authors state that L^1_loc-stability is verified but that uniqueness is not established.","tokens_in":40510,"tokens_out":5175,"duration_ms":67775,"significance":"If the classification is correct, it is a substantial advance over the authors' previous results, which were restricted to right states near the oil vertex or inside the quadrilateral O-E-U-D. The paper provides a detailed, systematic wave-curve construction for a large part of the saturation triangle in a physically motivated viscosity regime, with explicit claims for each subregion and supporting numerical simulations. Strengths include the use of an explicit Corey model that permits concrete computations, a transparent statement that uniqueness is not proved, and numerical validation in Section 5. However, the advertised scope is broader than what is actually established: the admissibility analysis is carried out for D(S)=I, the single-branch restriction on the double contact locus is not quantified, and the claimed L^1_loc-stability is not demonstrated in the text. These issues affect the central classification claim and need to be addressed before the paper can be accepted as a definitive classification.","major_comments":[{"comment":"The viscous profile admissibility criterion is implemented only for D(S)=I, as stated in Section 3.1: 'we assume D(S) to be the identity matrix and numerically verify the admissibility of shocks.' All admissible shock segments used in Claims 4.1-4.30, such as [A1,R), [A2,A3], and (R,A1], are selected by numerical checks of the D=I traveling-wave ODE. For a physical positive-definite capillary pressure matrix that is not proportional to the identity, the ODE vector field in Eq. (9) changes, the invariance of the segments G-D, W-E, and O-B is lost, and the criterion based on 'relative positions of M and N with respect to the corresponding segments' no longer applies. The admissible Hugoniot set, and hence the backward fast wave curve W_f(R), can change. Since the abstract and conclusion present the classification for three-phase flow in porous media without this qualification, the central claim is currently only supported for the special artificial case D=I. The authors should either prove that the admissible set is invariant under all positive-definite D in the allowed class, or explicitly restrict all theorems and claims to D=I and adjust the abstract and conclusion accordingly.","section":"Section 3.1, Eq. (9), Remark 3.1"},{"comment":"The paper imposes an additional restriction on the viscosities, namely that the double contact locus has only one branch in the saturation triangle, but this restriction is never translated into explicit inequalities on r_w and r_g. The text says 'The viscosities are further restricted so that the double contact locus, see Definition 3.5, possesses only one branch' and that the boundaries AO and A'O are defined by this restriction, yet no formula or proof is given. The numerical computations use the single parameter set (µ_w=1, µ_o=9.5, µ_g=0.45). Consequently the statement in the abstract that the classification 'remains valid for all viscosity variations satisfying the inequalities (5)' is not established, because the single-branch condition is an additional, unquantified hypothesis that can fail even when (5) holds. The authors need to provide the explicit parameter region or reformulate the main theorem to cover only the cases for which the topological assumptions are verified.","section":"Section 2, Eq. (5)-(6), Definition 3.5"},{"comment":"The paper claims L^1_loc-stability of the Riemann solution with respect to variations in the data, but no definition, theorem, proof, or numerical stability study is provided. Section 5 contains three comparisons of analytical and numerical saturation profiles at t_D=1; these single-time comparisons do not constitute a verification of L^1_loc-stability, which concerns continuous dependence of the solution as a curve in L^1_loc on the initial data. Either a precise stability statement with supporting argument or numerical experiments measuring data-to-solution continuity should be added, or the claim of L^1_loc-stability should be removed from the abstract, introduction, and conclusion.","section":"Abstract, Conclusion, Section 5"},{"comment":"The central classification is presented as a series of claims whose support is largely numerical or visual. For instance, Claim 4.1 and subsequent claims describe the complete backward fast wave curve W_f(R) by referring to figures and rely on statements such as 'We verify numerically that ... satisfy the viscous profile admissibility criterion.' The orderings of the intersection points L_1, L_R, L_*, L_3, L_2 along the edge G-W, which determine the left-state intervals in each claim, are asserted from the figures rather than proved. Combined with the explicit statement that uniqueness is not established, the phrase 'we classify all Riemann solution problems' overstates what is demonstrated: the paper provides a conjectured classification based on numerical evidence and geometric inspection. The authors should either supply proofs for the structural properties of the wave curves and the left-state orderings, or clearly rephrase the claims as a numerically supported classification and remove the word 'all' from the main claims.","section":"Section 4, Claims 4.1-4.30"}],"minor_comments":[{"comment":"The claim states that R lies in subregion Θ^b_1, but the surrounding text and Figure 13 indicate that the intended subregion is Θ^e_1; this should be corrected.","section":"Claim 4.5"},{"comment":"The caption refers to 'the ODE system (4.3)', but the relevant system is Eq. (9); the reference should be updated.","section":"Figure 25 caption"},{"comment":"There are typos such as 'ans-rarefaction' instead of 'an s-rarefaction' and 'The numerical simulation match' instead of 'The numerical simulation matches'; these should be corrected.","section":"Section 5, Cases 1 and 3"},{"comment":"The definition of ξ contains an unbalanced parenthesis: 'ξ=x−σt)/ε' should be 'ξ=(x−σt)/ε'.","section":"Section 3.1, Eq. (9) and surrounding text"},{"comment":"The sentence 'Refer to Fig. 6W f (R) comprises states M...' is missing a period after 'Fig. 6'; it should read 'Refer to Fig. 6. W_f(R) comprises states M...'.","section":"Claim 4.1"}],"recommendation":"major_revision","confidential_remarks":"The paper is a serious and detailed extension of the authors' program, and the overall construction appears internally consistent for the D=I model with the specific parameter values used. My main concern is that the advertised scope—classification for three-phase flow in porous media, valid for all viscosity variations satisfying (5), with L^1_loc-stability—exceeds what is actually proved. If the authors are willing to restrict the main theorems to D=I and to the explicitly verified parameter regime, and to either prove or clearly label the stability and uniqueness limitations, the paper could become acceptable. Given the load-bearing nature of the D=I assumption and the unquantified double-contact restriction, I recommend major revision rather than rejection."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"The real news here is the scope: previous work covered the corner O or the quadrilateral O-E-U-D, and this paper pushes the right-state classification to nearly the whole saturation triangle for the heavy-oil viscosity regime. That is a genuine extension, and the amount of wave-curve geometry involved is considerable. Thirty claims, each with explicit wave sequences and boundary states, plus three numerical simulations that match the analytical profiles well. If you work on hyperbolic models of multiphase flow, this is a useful reference to have on hand.\n\nThe paper is also unusually honest about its own limits. It states plainly that uniqueness is not proved, that L1_loc-stability is verified numerically rather than established analytically, and that shock admissibility is checked with the viscosity matrix taken to be the identity. That last point is the one to pay attention to. It appears in Section 3.1 and Remark 3.1, and it is load-bearing: the admissible shock segments, and therefore the wave curves and the classification itself, depend on the D=I ODE. The abstract and conclusion do not mention this simplification, so the advertised claim \"three-phase flow in porous media\" overreaches. For a non-identity positive-definite capillary pressure matrix, the invariant segments G-D, W-E, O-B are no longer invariant, and the shock admissibility criterion changes. That does not invalidate the paper as a classification of the D=I model, but it does mean the physical interpretation is narrower than the abstract suggests.\n\nTwo smaller soft spots. First, the L1_loc-stability statement is a numerical check, not a proof; the authors say so, but the abstract's \"We verify\" may read as stronger. Second, the claims are supported by figures and simulations rather than by rigorous proofs of admissibility for every shock segment. That is normal practice for this kind of wave-curve analysis, and the internal consistency of the construction is credible, but it is not formal verification.\n\nBottom line: this deserves a serious referee. The classification is a real step forward within an established research program, and the limitations are stated in the body even if the abstract glosses over them. I would send it to review with the request that the authors either justify the D=I assumption physically or qualify the claims in the abstract and conclusion. I would also ask them to separate \"verified numerically\" from \"proved\" in the stability statement. Not a desk reject; a revise-and-resubmit with clear requests.","headline":"A substantial extension of the three-phase Riemann classification, honestly limited by the D=I admissibility assumption and by numerical rather than proof-based verification.","tokens_in":41001,"tokens_out":1387,"would_cite":true,"duration_ms":25947,"reading_group":"maybe","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["35L65","76S05","76T30"],"pacs":[],"model":"deepseek-v4-flash","headline":"For heavy-oil reservoirs, the paper classifies every Riemann solution connecting a water–gas injection state on the edge G–W to almost any displaced state in the saturation triangle, and shows the resulting wave pattern is stable to small…","keywords":["Riemann solutions","multiphase flow in porous media","heavy oil","wave curve method","viscous profile admissibility","umbilic point","saturation triangle","water-alternating-gas injection"],"falsifier":"Pick a viscosity ratio satisfying the inequalities in (5), replace D(S) = I by a positive-definite anisotropic matrix in the parabolic system (8), and recompute the viscous-profile condition for the shock segment [A_2, A_3] that is declared admissible for R in region Θ_1; if that segment loses its profile or a previously non-admissible segment gains one, the classification is specific to the identity-matrix model rather than to the physical three-phase flow.","tokens_in":40076,"feed_emoji":"🛢️","tokens_out":8570,"duration_ms":91665,"temperature":0.7,"pith_summary":"This paper treats the injection of water and gas into a heavy-oil reservoir as a Riemann problem for a system of two conservation laws, and tries to determine, for every initial mixture of oil, water, and gas, exactly which sequence of waves—rarefactions, shocks, and composite waves—will form. The main claim is a full classification for the case where the injected state lies on the water–gas edge of the saturation triangle and the displaced state lies anywhere in the triangle outside two small boundary regions. The classification holds for viscosity ratios satisfying two explicit inequalities, precisely the regime in which the umbilic point sits close to the oil vertex, and it uses only classical waves. The paper also proves that the constructed solutions depend continuously on the initial data in the locally integrable sense, and it corroborates the predicted profiles with direct numerical simulation. If the classification is right, it gives a predictive catalog for water-alternating-gas recovery in heavy oil and a template for more general permeability models.","feed_headline":"All Riemann solutions mapped for heavy-oil three-phase flow","feed_subtitle":"Injection states on the gas-water edge now cover nearly the whole saturation triangle, with stability verified.","key_machinery":"The wave curve method. For each right state R, one constructs the backward fast wave curve W_f(R) out of fast rarefaction segments, admissible fast shock segments of the Hugoniot locus, and composite segments defined by extensions of rarefaction and inflection loci; for each left state L on G–W, one constructs the forward slow wave curve and requires speed compatibility (end of slow group no faster than start of fast group) to assemble the solution. The organization of the classification is carried out by bifurcation loci in the saturation triangle—the secondary bifurcation loci E–W, G–D, O–B, the inflection loci, the hysteresis loci, the double-contact and mixed-contact loci, and the tangential extension locus T_I—whose intersections define regions of right states with a fixed solution structure.","core_discovery":"The central discovery is that in the heavy-oil viscosity regime, every Riemann solution with left state on the edge G–W and right state in almost all of the saturation triangle consists of at most two wave groups—a slow group followed by a fast group—with no undercompressive or overcompressive shocks. The paper partitions the saturation triangle into regions Λ, Θ, Ω, Γ, subdivided by bifurcation loci (inflection, hysteresis, double-contact, mixed-contact, and tangential extensions), and for each subregion it lists the admissible backward fast wave curve and the resulting composition path for every left state L on G–W. Shock admissibility is decided by the viscous profile criterion under the simplification that the viscosity matrix is the identity, and each admissible segment is verified numerically. The paper verifies L1_loc stability of the Riemann solution with respect to variations in the data and presents numerical simulations that match the analytical saturation profiles.","pith_inferences":["A natural next step is the two small unclassified strips near the G–O and W–O edges; the detached Hugoniot branches that are non-admissible under the identity viscosity matrix may become admissible under other capillary matrices, so those strips could harbor undercompressive or overcompressive waves even in this viscosity regime.","The identity-matrix simplification is probably the most fragile link: a diagonal but non-scalar capillary matrix, let alone an anisotropic one, would change the traveling-wave ODEs and could eliminate or create admissible shock segments, so the classification should be rechecked numerically for at least one physically motivated non-identity D(S).","The authors' restriction to one branch of the double-contact locus (Definition 3.5) suggests that for viscosity ratios satisfying (5) but closer to equality, additional contact branches enter the triangle and the subdivision of Θ, Ω, Γ would need to be refined; this gives a concrete parameter-space boundary where the present classification breaks.","The stability proved is L1_loc with respect to data variations within the class; if a full uniqueness theorem were desired, one would need to compare against an entropy or front-tracking selection rule, and a numerical search for alternative paths exactly at subregion boundaries (e.g., L = L* cases) could test whether the triple-shock rule really collapses the two candidate paths to the same solut"],"forward_implications":["In the heavy-oil regime, water-alternating-gas injection into almost any oil-water-gas mixture produces only classical shock and rarefaction waves, so front-tracking simulators can use the explicit composition paths instead of costly Riemann solvers.","The classification provides a complete map of which right states reach the oil vertex region versus stay near the gas-water edge, which bears on oil displacement efficiency for given injection compositions.","The L1_loc stability result means that small perturbations in measured initial saturations lead to small perturbations of the predicted saturation profiles, supporting the use of these solutions for uncertainty quantification.","The framework reduces the previously treated regions (near oil vertex and quadrilateral O–E–U–D) to a nearly global picture of the saturation triangle, leaving only two small boundary strips unclassified.","Because the model uses quadratic Corey permeabilities, the same wave-curve construction can be repeated for other Corey exponents, so the paper supplies the base case for a wider family of heavy-oil relative-permeability laws."],"supporting_citations":[{"why":"Establishes the earlier Riemann classification for left states on G–W with right states near the oil vertex; the present work extends that result.","marker":"[1]"},{"why":"Proved that only classical waves occur in the quadrilateral O–E–U–D under inequalities (5); the present work extends this to nearly the whole triangle.","marker":"[2]"},{"why":"Introduced the wave curve method solution for three-phase flow in virgin reservoirs, the method used throughout this paper.","marker":"[3]"},{"why":"Establishes uniqueness and properties of characteristic speeds (real, positive, with the umbilic point), a foundation for the wave-curve construction.","marker":"[4]"},{"why":"Interactive graphical Riemann problem solver used for numerical verification of admissible shock segments and of the predicted solutions.","marker":"[10]"},{"why":"D.Sc. thesis establishing diffusive effects and the identity-viscosity-matrix simplification used for shock admissibility.","marker":"[11]"},{"why":"Generalized the wave curve method to singular Riemann problems; this generalization underpins the construction of composite and extension loci.","marker":"[13]"},{"why":"The original wave curve method for 2x2 conservation laws, which the paper adapts.","marker":"[18]"},{"why":"Schaeffer-Shearer classification of umbilic points; the viscosity regime considered corresponds to Case II, locating the umbilic point close to the oil vertex.","marker":"[28]"}],"fun_headline_variants":["Heavy-oil flow: every Riemann solution classified","Riemann map covers nearly all saturation states for heavy oil","Stable Riemann solutions for heavy-oil gas-water injection","Full Riemann solution classification for heavy-oil displacement","Three-phase flow Riemann solutions: heavy-oil case fully solved"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The load-bearing premise is that shock admissibility can be decided by viscous profiles with the capillary-pressure matrix D(S) taken to be the identity matrix; if the true capillary matrix is not a scalar multiple of the identity, the set of admissible shocks—and hence the wave sequences in the classification—could change.","fun_headline_variants_meta":{"raw":{"variants":["Heavy-oil flow: every Riemann solution classified","Riemann map covers nearly all saturation states for heavy oil","Stable Riemann solutions for heavy-oil gas-water injection","Full Riemann solution classification for heavy-oil displacement","Three-phase flow Riemann solutions: heavy-oil case fully solved"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000996,"raw_usage":{"total_tokens":4216,"prompt_tokens":941,"completion_tokens":3275,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":557,"completion_tokens_details":{"reasoning_tokens":3196}},"tokens_in":557,"tokens_out":3275,"duration_ms":27514,"temperature":1.0,"reasoning_tokens":3196,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-07T05:12:09.138178+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Pick a viscosity ratio satisfying the inequalities in (5), replace D(S) = I by a positive-definite anisotropic matrix in the parabolic system (8), and recompute the viscous-profile condition for the shock segment [A_2, A_3] that is declared admissible for R in region Θ_1; if that segment loses its profile or a previously non-admissible segment gains one, the classification is specific to the identity-matrix model rather than to the physical three-phase flow.","supporting_citations":[{"cited_title":"Andrade, A","cited_arxiv_id":null,"evidence_quote":"Establishes the earlier Riemann classification for left states on G–W with right states near the oil vertex; the present work extends that result."},{"cited_title":"Andrade, A","cited_arxiv_id":null,"evidence_quote":"Proved that only classical waves occur in the quadrilateral O–E–U–D under inequalities (5); the present work extends this to nearly the whole triangle."},{"cited_title":"Azevedo, A","cited_arxiv_id":null,"evidence_quote":"Introduced the wave curve method solution for three-phase flow in virgin reservoirs, the method used throughout this paper."},{"cited_title":"Azevedo, A","cited_arxiv_id":null,"evidence_quote":"Establishes uniqueness and properties of characteristic speeds (real, positive, with the umbilic point), a foundation for the wave-curve construction."},{"cited_title":"ELI , I nteractive G raphical R iemann P roblem S olver","cited_arxiv_id":null,"evidence_quote":"Interactive graphical Riemann problem solver used for numerical verification of admissible shock segments and of the predicted solutions."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"D.Sc. thesis establishing diffusive effects and the identity-viscosity-matrix simplification used for shock admissibility."},{"cited_title":"Isaacson, D","cited_arxiv_id":null,"evidence_quote":"Generalized the wave curve method to singular Riemann problems; this generalization underpins the construction of composite and extension loci."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"The original wave curve method for 2x2 conservation laws, which the paper adapts."},{"cited_title":"The classification of 2 2 systems of non-strictly hyperbolic conservation laws, with application to oil recovery","cited_arxiv_id":null,"evidence_quote":"Schaeffer-Shearer classification of umbilic points; the viscosity regime considered corresponds to Case II, locating the umbilic point close to the oil vertex."}],"review_version":1}