ERROR neqsim.thermo.phase.PhaseEos - Failed to solve for molarVolume within the iteration limit.

[sface-admin]

I often get this error where composition has pseudo components. Is this particularly an error due to python envelope?
The snip of example code below.

from neqsim.thermo.thermoTools import *
import matplotlib.pyplot as plt

fluid1 = fluid("srk") # create a fluid using the SRK-EoS
fluid1.setTemperature(80, "C")
fluid1.setPressure(50, "bara")
fluid1.addComponent("nitrogen", 0)
fluid1.addComponent("CO2", 1.005)
fluid1.addComponent("H2S", 0)
fluid1.addComponent("methane", 24.62)
fluid1.addComponent("ethane", 2.704)
fluid1.addComponent("propane", 2.706)
fluid1.addComponent("i-butane", 1.309)
fluid1.addComponent("n-butane", 2.031)
fluid1.addComponent("i-pentane", 2.047)
fluid1.addComponent("n-pentane", 1.49)
fluid1.addTBPfraction("C6", 3.176, 85.3/1000.0, 0.6649)
fluid1.addTBPfraction("C7", 5.331, 93.2/1000.0, 0.7266)
fluid1.addTBPfraction("C8", 7.796, 105.1/1000.0, 0.7607)
fluid1.addTBPfraction("C9", 5.179, 119.4/1000.0, 0.7812)
fluid1.addPlusFraction("C10+", 40.616, 268/1000.0, 0.865)
fluid1.setMixingRule("classic")
fluid1.setMultiPhaseCheck(True)
fluid1.useVolumeCorrection(True)
TPflash(fluid1)
fluid1.getInterphaseProperties().setInterfacialTensionModel("gas", "oil", "Parachor")
fluid1.initProperties()
fluid1.display()
thermoOps = jNeqSim.thermodynamicOperations.ThermodynamicOperations(fluid1)
thermoOps.calcPTphaseEnvelope()
plt.plot(list(thermoOps.getOperation().get("dewT2")), list(thermoOps.getOperation().get("dewP2")), label="DewPoint")
plt.plot(list(thermoOps.getOperation().get("bubT")), list(thermoOps.getOperation().get("bubP")), label="BubblePoint")
plt.show()

[EvenSol]

Thanks for pointing at this.

It is a weakness in the phase envelope algorithm not accepting components with mole fraction 0. We will fix this in a later release, and mark this as a bug.

To run the script without any error, you can skip adding components where number of moles are 0.
The script could be:

from neqsim.thermo.thermoTools import *
import matplotlib.pyplot as plt

fluid1 = fluid("srk") # create a fluid using the SRK-EoS
fluid1.setTemperature(80, "C")
fluid1.setPressure(50, "bara")
fluid1.addComponent("CO2", 1.005)
fluid1.addComponent("methane", 24.62)
fluid1.addComponent("ethane", 2.704)
fluid1.addComponent("propane", 2.706)
fluid1.addComponent("i-butane", 1.309)
fluid1.addComponent("n-butane", 2.031)
fluid1.addComponent("i-pentane", 2.047)
fluid1.addComponent("n-pentane", 1.49)
fluid1.addTBPfraction("C6", 3.176, 85.3/1000.0, 0.6649)
fluid1.addTBPfraction("C7", 5.331, 93.2/1000.0, 0.7266)
fluid1.addTBPfraction("C8", 7.796, 105.1/1000.0, 0.7607)
fluid1.addTBPfraction("C9", 5.179, 119.4/1000.0, 0.7812)
fluid1.addPlusFraction("C10+", 40.616, 268/1000.0, 0.865)
fluid1.setMixingRule("classic")
fluid1.setMultiPhaseCheck(True)
fluid1.useVolumeCorrection(True)
TPflash(fluid1)
fluid1.getInterphaseProperties().setInterfacialTensionModel("gas", "oil", "Parachor")
fluid1.initProperties()
fluid1.display()
thermoOps = jNeqSim.thermodynamicOperations.ThermodynamicOperations(fluid1)
thermoOps.calcPTphaseEnvelope()
plt.plot(list(thermoOps.getOperation().get("dewT")), list(thermoOps.getOperation().get("dewP")), label="DewPoint")
plt.plot(list(thermoOps.getOperation().get("bubT")), list(thermoOps.getOperation().get("bubP")), label="BubblePoint")
plt.show()

[sface-admin]

Thanks, that does solve for this composition. However I get similar error for another gas composition, where all components are non zero. Notice that error starts with addition of C7 pseudo onwards. See the script below,

fluid1 = fluid("srk")
fluid1.setTemperature(40, "C")
fluid1.setPressure(50, "bara")
fluid1.addComponent("nitrogen", 0.88)
fluid1.addComponent("CO2", 5.7)
fluid1.addComponent("methane", 86.89)
fluid1.addComponent("ethane", 3.59)
fluid1.addComponent("propane", 1.25)
fluid1.addComponent("i-butane", 0.19)
fluid1.addComponent("n-butane", 0.35)
fluid1.addComponent("i-pentane", 0.12)
fluid1.addComponent("n-pentane", 0.12)
fluid1.addTBPfraction("C6", 0.15, 86/1000.0, 0.672)
fluid1.addTBPfraction("C7", 0.2, 96/1000.0, 0.737)
fluid1.addTBPfraction("C8", 0.22, 106/1000.0, 0.767)
fluid1.addTBPfraction("C9", 0.13, 121/1000.0, 0.783)
fluid1.addPlusFraction("C10+", 0.21, 172/1000.0, 0.818)
fluid1.setMixingRule("classic")
fluid1.setMultiPhaseCheck(True)
fluid1.useVolumeCorrection(True)
initPhysicalProperties()
thermoOps = jNeqSim.thermodynamicOperations.ThermodynamicOperations(fluid1)
thermoOps.calcPTphaseEnvelope()
plt.plot(list(thermoOps.getOperation().get("dewT")), list(thermoOps.getOperation().get("dewP")), label="DewPoint")
plt.plot(list(thermoOps.getOperation().get("bubT")), list(thermoOps.getOperation().get("bubP")), label="BubblePoint")
plt.show()

[EvenSol]

All this seems to work in the new version and using the method:

thermoOps = jNeqSim.thermodynamicOperations.ThermodynamicOperations(fluid1)
thermoOps.calcPTphaseEnvelope2()

See the tests in the notebook:
https://colab.research.google.com/drive/1tuSXN5Myk97SSxxJJ21vo4G8lwILT9K2#scrollTo=iKZ_5EKuf09H