Class PhaseSrkCPAreduced

java.lang.Object
neqsim.thermo.phase.Phase
neqsim.thermo.phase.PhaseEos
neqsim.thermo.phase.PhaseSrkEos
neqsim.thermo.phase.PhaseSrkCPA
neqsim.thermo.phase.PhaseSrkCPAs
neqsim.thermo.phase.PhaseSrkCPAreduced
All Implemented Interfaces:
Serializable, Cloneable, PhaseCPAInterface, PhaseEosInterface, PhaseInterface, ThermodynamicConstantsInterface
public class PhaseSrkCPAreduced extends PhaseSrkCPAs
CPA phase class with site symmetry reduction and Broyden quasi-Newton acceleration.

Reduces the coupled (n_s + 1)-dimensional Newton system to (p + 1) dimensions, where p is the number of unique association site types. Sites on the same component with identical bonding patterns (same delta row in the association matrix) are grouped into a single "type" with a multiplicity factor. This exploits the mathematical theorem that equivalent sites converge to equal site fractions at equilibrium.

Combined with Broyden rank-1 updates of the inverse Jacobian after initial convergence, this yields compounded speedup from both dimension reduction and reduced per-iteration cost.

Dimension reduction examples:

Dimension reduction for common CPA systems
System n_s (full) p (reduced) Jacobian speedup
Pure water (4C) 4 2 4.6x
Water + methanol (4C+2B) 6 4 2.7x
Water + MEG (4C+4C) 8 4 5.8x
NG + water + TEG 8 4 5.8x
Version:
1.0
Author:
Even Solbraa
See Also:

  • Field Details

    • serialVersionUID

      private static final long serialVersionUID
      Serialization version UID.
      See Also:

    • MAX_ITERATIONS

      static final int MAX_ITERATIONS
      Maximum iterations for the coupled Newton/Broyden solver.
      See Also:

    • CONVERGENCE_TOL

      static final double CONVERGENCE_TOL
      Convergence tolerance for the max residual.
      See Also:

    • MIN_NEWTON_STEPS

      static final int MIN_NEWTON_STEPS
      Minimum full Newton steps before switching to Broyden.
      See Also:

    • BROYDEN_SWITCH_TOL

      static final double BROYDEN_SWITCH_TOL
      Max residual below which Broyden updates are activated.
      See Also:

    • STALL_RATIO

      static final double STALL_RATIO
      Ratio threshold for detecting stalling.
      See Also:

    • MAX_REL_STEP

      static final double MAX_REL_STEP
      Maximum relative step.
      See Also:

    • numTypes

      private transient int numTypes
      Number of unique site types for the current system.

    • siteToType

      private transient int[] siteToType
      Maps individual site index to type index.

    • typeRepSite

      private transient int[] typeRepSite
      Representative individual site index for each type.

    • typeMult

      private transient int[] typeMult
      Multiplicity (count of equivalent sites) for each type.

    • typeCompIdx

      private transient int[] typeCompIdx
      Component index for each type.

    • callCount

      long callCount
      Number of molarVolume calls.

    • totalIters

      long totalIters
      Total iterations across all calls.

    • jacobianEvals

      long jacobianEvals
      Number of full Jacobian evaluations.

    • broydenUpdates

      long broydenUpdates
      Number of Broyden rank-1 updates.

    • fallbackCount

      long fallbackCount
      Number of solver fallbacks to the base class.

    • lastNumTypes

      int lastNumTypes
      Number of site types (last call).

    • lastFullSites

      int lastFullSites
      Number of full sites (last call).

    • workKlk

      private transient double[][] workKlk

    • workHess

      private transient double[][] workHess

    • workKsi

      private transient double[] workKsi

    • workM

      private transient double[] workM

    • workKlkKsi

      private transient double[] workKlkKsi

    • workXV

      private transient double[] workXV

    • workP

      private transient int workP

  • Constructor Details

    • PhaseSrkCPAreduced

      public PhaseSrkCPAreduced()
      Construct a PhaseSrkCPAreduced phase.

  • Method Details

    • clone

      public PhaseSrkCPAreduced clone()

      clone.

      Specified by:
      clone in interface PhaseInterface
      Overrides:
      clone in class PhaseSrkCPAs
      Returns:
      a PhaseInterface object

    • addComponent

      public void addComponent(String name, double moles, double molesInPhase, int compNumber)

      Add component to component array and update moles variables.

      Specified by:
      addComponent in interface PhaseInterface
      Overrides:
      addComponent in class PhaseSrkCPAs
      Parameters:
      name - Name of component.
      moles - Total number of moles of component.
      molesInPhase - Number of moles in phase.
      compNumber - Index number of component in phase object component array.

    • molarVolume

      public double molarVolume(double pressure, double temperature, double A, double B, PhaseType pt) throws IsNaNException, TooManyIterationsException

      molarVolume.

      Reduced-dimension coupled Newton/Broyden molar volume solver. Exploits association site symmetry to work in (p+1) dimensions where p is the number of unique site types, combined with Broyden rank-1 inverse-Jacobian updates after the initial Newton phase.

      Specified by:
      molarVolume in interface PhaseInterface
      Overrides:
      molarVolume in class PhaseSrkCPA
      Parameters:
      pressure - a double
      temperature - a double
      A - a double
      B - a double
      pt - the PhaseType of the phase
      Returns:
      a double
      Throws:
      IsNaNException - if any.
      TooManyIterationsException - if any.

    • buildSiteTypeMap

      private void buildSiteTypeMap(int ns)
      Build the site type map by grouping equivalent association sites.

      Two individual sites are equivalent if they belong to the same component and have identical deltaNog rows (same bonding pattern to all other sites). This corresponds to sites with the same charge in the CPA association scheme — e.g., the two electron-donor sites of water (4C) have identical interactions with all other sites.

      Parameters:
      ns - total number of individual association sites

    • expandAndSetSiteFractions

      private void expandAndSetSiteFractions(double[] xType, int ns)
      Expand reduced site fraction values to all individual sites on components.
      Parameters:
      xType - reduced site fraction array (length p)
      ns - total number of individual sites

    • readXsiteFromComponents

      private void readXsiteFromComponents(double[] xSite, int ns)
      Read individual site fractions from component objects into a flat array.
      Parameters:
      xSite - array to fill (length ns)
      ns - total number of individual sites

    • buildReducedJacobian

      private void buildReducedJacobian(double[][] jac, double[] xType, double[] tMoles, double[] redSum, double totalVol, double gdv1, double zeta, double btemp, int dim, int p)
      Build the analytic Jacobian for the reduced (p+1)-dimensional system.

      Block structure:

      • J_XX[a][b]: site fraction block with multiplicity weighting
      • j_X_zeta[a]: site fraction sensitivity to volume (via zeta)
      • j_zeta_X[a]: volume equation sensitivity to site fractions
      • j_zeta_zeta: volume equation self-sensitivity
      Parameters:
      jac - Jacobian matrix to populate (dim x dim)
      xType - reduced site fraction values (length p)
      tMoles - moles per type (length p)
      redSum - precomputed reduced summation for each type (length p)
      totalVol - total volume V
      gdv1 - g'(V) - 1/V
      zeta - current B/(nV)
      btemp - co-volume B
      dim - system dimension (p+1)
      p - number of unique site types

    • updateDeltaWithG

      private void updateDeltaWithG(int ns)
      Multiply delta by g-function for all individual sites.
      Parameters:
      ns - total number of individual sites

    • invertMatrix

      private static boolean invertMatrix(double[][] a, double[][] ainv, int n)
      Invert matrix A into Ainv using Gauss-Jordan elimination with partial pivoting.
      Parameters:
      a - input matrix (will be modified by the elimination)
      ainv - output inverse matrix
      n - matrix dimension
      Returns:
      true if successful, false if singular

    • broydenUpdate

      private static void broydenUpdate(double[][] invJ, double[] dxVec, double[] dfVec, int n)
      Broyden "good" rank-1 update of the inverse Jacobian H.

      H_new = H + (dx - H*df) * (dx^T * H) / (dx^T * H * df). Cost: O(n^2).

      Parameters:
      invJ - inverse Jacobian (updated in-place)
      dxVec - step vector x_{k+1} - x_k
      dfVec - residual change R(x_{k+1}) - R(x_k)
      n - system dimension

    • ensureWorkArrays

      private void ensureWorkArrays(int p)
      Ensure reduced work arrays are allocated.
      Parameters:
      p - number of unique site types

    • initCPAMatrix

      public void initCPAMatrix(int type)

      initCPAMatrix.

      Override type 1 initialization to use reduced-dimension site type computation. Higher-order volume derivatives (FCPA, dFCPAdV, dFCPAdVdV, dFCPAdVdVdV) are computed using the type grouping, reducing linear system size from n_s to p.

      Overrides:
      initCPAMatrix in class PhaseSrkCPA
      Parameters:
      type - a int

    • solveLinearSystem

      private static void solveLinearSystem(double[][] a, double[] b, int n)
      Solve the linear system A*x = b in-place (b overwritten with solution) using Gaussian elimination with partial pivoting.
      Parameters:
      a - coefficient matrix (modified in-place)
      b - right-hand side vector (overwritten with solution)
      n - system dimension

    • molarVolumeChangePhase

      public double molarVolumeChangePhase(double pressure, double temperature, double A, double B, PhaseType pt) throws IsNaNException, TooManyIterationsException

      molarVolumeChangePhase.

      Overrides:
      molarVolumeChangePhase in class PhaseSrkCPA
      Parameters:
      pressure - a double
      temperature - a double
      A - a double
      B - a double
      pt - the PhaseType of the phase
      Returns:
      a double
      Throws:
      IsNaNException - if any.
      TooManyIterationsException - if any.

    • getCallCount

      public long getCallCount()
      Get the total number of molarVolume calls.
      Returns:
      call count

    • getTotalIters

      public long getTotalIters()
      Get the total Newton/Broyden iterations across all calls.
      Returns:
      total iterations

    • getJacobianEvals

      public long getJacobianEvals()
      Get the total number of full Jacobian evaluations.
      Returns:
      Jacobian evaluation count

    • getBroydenUpdates

      public long getBroydenUpdates()
      Get the total number of Broyden rank-1 updates.
      Returns:
      Broyden update count

    • getFallbackCount

      public long getFallbackCount()
      Get the number of fallbacks to the base class solver.
      Returns:
      fallback count

    • getLastNumTypes

      public int getLastNumTypes()
      Get the number of unique site types from the last molarVolume call.
      Returns:
      number of unique types (p)

    • getLastFullSites

      public int getLastFullSites()
      Get the total number of individual sites from the last molarVolume call.
      Returns:
      number of full sites (n_s)

    • resetProfiling

      public void resetProfiling()
      Reset all profiling counters to zero.