Skip to content

Commit

Permalink
feat: update pipeBeggsAndBrills (#1273)
Browse files Browse the repository at this point in the history
  • Loading branch information
@Sviatose
Sviatose authored Feb 9, 2025
1 parent 155e8ea commit ec372df
Show file tree
Hide file tree
Showing 2 changed files with 255 additions and 6 deletions.
198 changes: 193 additions & 5 deletions src/main/java/neqsim/process/equipment/pipeline/PipeBeggsAndBrills.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -146,6 +146,40 @@ public class PipeBeggsAndBrills extends Pipeline {
private List<Double> elevationProfile;
private List<Integer> incrementsProfile;

// Flag to run isothermal calculations
private boolean runAdiabatic = true;
private boolean runConstantSurfaceTemperature = false;

private double constantSurfaceTemperature;

private double heatTransferCoefficient;

private String heatTransferCoefficientMethod = "Estimated";


// Heat transfer parameters
double Tmi; // medium temperature
double Tmo; // outlet temperature
double Ts; // wall temperature
double error; // error in heat transfer
double iterationT; // iteration in heat transfer
double dTlm; // log mean temperature difference
double cp; // heat capacity
double q1; // heat transfer
double q2;
double ReNoSlip; //
double S = 0;
double rhoNoSlip = 0;
double muNoSlip = 0;
double thermalConductivity;
double Pr; // Prandtl number
double frictionFactor;
double frictionTwoPhase;
double Nu;
double criticalPressure;
double hmax;
double X;

/**
* Constructor for PipeBeggsAndBrills.
*
Expand Down Expand Up @@ -302,6 +336,40 @@ public void setRunIsothermal(boolean runIsothermal) {
this.runIsothermal = runIsothermal;
}

/**
* <p>
* Setter for the field <code>constantSurfaceTemperature</code>.
* </p>
*
* @param temperature a double
*/
public void setConstantSurfaceTemperature(double temperature, String unit) {
if (unit.equals("K")) {
this.constantSurfaceTemperature = temperature;
} else if (unit.equals("C")) {
this.constantSurfaceTemperature = temperature + 273.15;
} else {
throw new RuntimeException("unit not supported " + unit);
}
this.runIsothermal = false;
this.runAdiabatic = false;
this.runConstantSurfaceTemperature = true;
}


/**
* <p>
* Setter for the field <code>heatTransferCoefficient</code>.
* </p>
*
* @param heatTransferCoefficient a double
*/
public void setHeatTransferCoefficient(double heatTransferCoefficient) {
this.heatTransferCoefficient = heatTransferCoefficient;
this.heatTransferCoefficientMethod = "Defined";
}


/**
* Converts the input values from the system measurement units to imperial units. Needed because
* the main equations and coefficients are developed for imperial system
Expand Down Expand Up @@ -650,16 +718,16 @@ public double calcFrictionPressureLoss() {
mixtureViscosityProfile.add(muNoSlip);
mixtureDensityProfile.add(rhoNoSlip * 16.01846);

double ReNoSlip =
ReNoSlip =
rhoNoSlip * supMixVel * insideDiameter * (16 / (3.28 * 3.28)) / (0.001 * muNoSlip);

mixtureReynoldsNumber.add(ReNoSlip);

double E = pipeWallRoughness / insideDiameter;

// Haaland equation
double frictionFactor = Math.pow(1 / (-1.8 * Math.log10((E / 3.7) + (6.9 / ReNoSlip))), 2);
double frictionTwoPhase = frictionFactor * Math.exp(S);
frictionFactor = Math.pow(1 / (-1.8 * Math.log10((E / 3.7) + (6.9 / ReNoSlip))), 2);
frictionTwoPhase = frictionFactor * Math.exp(S);

frictionPressureLoss =
frictionTwoPhase * Math.pow(supMixVel, 2) * rhoNoSlip * (length) / (2 * insideDiameter);
Expand Down Expand Up @@ -741,6 +809,7 @@ public void run(UUID id) {
pressureDropProfile.add(pressureDrop);
pressureOut = inletPressure - pressureDrop;
pressureProfile.add(pressureOut);

if (pressureOut < 0) {
throw new RuntimeException(new neqsim.util.exception.InvalidOutputException(
"PipeBeggsAndBrills", "run: calcOutletPressure", "pressure out",
Expand All @@ -749,23 +818,133 @@ public void run(UUID id) {

system.setPressure(pressureOut);
if (!runIsothermal) {
testOps.PHflash(enthalpyInlet);
enthalpyInlet = calcHeatBalance(enthalpyInlet, system, testOps);
//testOps.PHflash(enthalpyInlet);
temperatureProfile.add(system.getTemperature());
} else {
testOps.TPflash();
}
system.initProperties();
temperatureProfile.add(system.getTemperature());
}
totalPressureDrop = pipeInletPressure - system.getPressure();
calcPressureDrop(); // to initialize final parameters
lengthProfile.add(cumulativeLength);
elevationProfile.add(cumulativeElevation);
incrementsProfile.add(getNumberOfIncrements());


outStream.setThermoSystem(system);
outStream.setCalculationIdentifier(id);
}

/**
* Estimates the heat transfer coefficient for the given system.
*
* @param system the thermodynamic system for which the heat transfer coefficient is to be estimated
* @return the estimated heat transfer coefficient
*/
public double estimateHeatTransferCoefficent(SystemInterface system){
cp = system.getCp("J/kgK");
thermalConductivity = system.getThermalConductivity();
Pr = 0.001*system.getViscosity("cP") * cp / thermalConductivity;
if( ReNoSlip < 3000 ){
Nu = 3.66;
}else{
if(Pr < 2000 && Pr > 0.5 && ReNoSlip < 5E6){
Nu = ((frictionTwoPhase/8)*(ReNoSlip - 1000)*Pr)/(1 + 12.7*Math.pow(frictionTwoPhase, 0.5)*(Math.pow(Pr,0.66)-1));
}
}
heatTransferCoefficient = Nu*thermalConductivity/(insideDiameter);


if (system.getNumberOfPhases() > 1){
X = system.getPhase(0).getFlowRate("kg/sec") / system.getFlowRate("kg/sec");
heatTransferCoefficient = (-31.469*Math.pow(X, 2) + 31.469*X+0.007)*heatTransferCoefficient;
// double Xtt = (Math.pow(((1-X)/X),0.9)*Math.pow((system.getPhase(0).getDensity()/mixtureLiquidDensity), 0.5)*
// Math.pow(mixtureLiquidViscosity/(0.001*system.getPhase(0).getViscosity()), 0.1));
// double E = 0.8897*(1/Xtt) + 1.0392;
// double Retp = ReNoSlip*Math.pow(E, 1.25);
// double S;
// if (Retp < 1E5){
// S = 0.9;
// }
// else if(Retp <= 1E6){
// S = -5.6*1E-7*Retp + 0.9405;
// }
// else if(Retp < 1E7)
// S = 1368.9*Math.pow(Retp, -0.601);
// else{
// S = 0.1;
// }

// double reducedPressure = system.getPressure()*100/criticalPressure;
// double Fp = 1.8*Math.pow(reducedPressure, 0.17) + 4*Math.pow(reducedPressure, 1.2) + 10*Math.pow(reducedPressure, 10);
// double hb = 0.00417*Math.pow(criticalPressure, 0.69)*Math.pow(heatFlux/100, 0.7)*Fp;
// heatTransferCoefficient = E*heatTransferCoefficient + S*hb;
}

return heatTransferCoefficient;
}

/**
* Calculates the temperature difference between the outlet and inlet of the system.
*
* @param system the thermodynamic system for which the temperature difference is to be calculated
* @return the temperature difference between the outlet and inlet
*/
public double calcTemperatureDifference(SystemInterface system){
double cp = system.getCp("J/kgK");
double Tmi = system.getTemperature("C");
double Ts = constantSurfaceTemperature - 273.15;
double TmoLower = Tmi; // Lower bound for Tmo
double TmoUpper = Ts; // Upper bound for Tmo
double Tmo = (TmoLower + TmoUpper) / 2; // Initial guess
double error = 999;
double tolerance = 0.01; // Tolerance for convergence
int maxIterations = 100; // Maximum number of iterations

if (heatTransferCoefficientMethod.equals("Estimated")) {
heatTransferCoefficient = estimateHeatTransferCoefficent(system);
}

for (int i = 0; i < maxIterations; i++) {
double dTlm = ((Ts - Tmo) - (Ts - Tmi)) / (Math.log((Ts - Tmo) / (Ts - Tmi)));
error = heatTransferCoefficient - system.getFlowRate("kg/sec") * cp * (Tmo - Tmi) / (3.1415 * insideDiameter * length * dTlm);

if (Math.abs(error) < tolerance) {
break; // Converged
}

if (error > 0) {
TmoLower = Tmo; // Adjust lower bound
} else {
TmoUpper = Tmo; // Adjust upper bound
}

Tmo = (TmoLower + TmoUpper) / 2; // New guess
}
return Tmo - Tmi;
}



/**
* Calculates the heat balance for the given system.
*
* @param enthalpy the initial enthalpy of the system
* @param system the thermodynamic system for which the heat balance is to be calculated
* @param testOps the thermodynamic operations to be performed
* @return the calculated enthalpy after performing the heat balance
*/
public double calcHeatBalance(double enthalpy, SystemInterface system, ThermodynamicOperations testOps){
double Cp = system.getCp("J/kgK");
if (!runAdiabatic) {
enthalpy = enthalpy + system.getFlowRate("kg/sec")*Cp*calcTemperatureDifference(system);
}
testOps.PHflash(enthalpy);
return enthalpy;
}

/**
* {@inheritDoc}
*
Expand Down Expand Up @@ -799,6 +978,15 @@ public double getInletSuperficialVelocity() {
/ (Math.PI / 4.0 * Math.pow(insideDiameter, 2.0));
}

/**
* Getter for the field <code>heatTransferCoefficient</code>.
*
* @return the heat transfer coefficient
*/
public double getHeatTransferCoefficient() {
return heatTransferCoefficient;
}

/**
* <p>
* getOutletSuperficialVelocity.
Expand Down
63 changes: 62 additions & 1 deletion src/test/java/neqsim/process/equipment/pipeline/BeggsAndBrillsPipeTest.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,7 @@ public void testFlowNoVolumeCorrection() {

ThermodynamicOperations testOps = new ThermodynamicOperations(testSystem);
testOps.TPflash();
testSystem.initPhysicalProperties();
// testSystem.initPhysicalProperties();

Assertions.assertEquals(testSystem.getPhase("oil").getFlowRate("m3/hr"),
testSystem.getFlowRate("m3/hr"), 1);
Expand Down Expand Up @@ -344,4 +344,65 @@ public void testPipeLineBeggsAndBrills4() {
Assertions.assertEquals(temperatureOut3, -8.81009355441591, 1);
Assertions.assertEquals(pressureOut3, 18.3429, 1);
}


@Test
public void testPipeLineBeggsAndBrills5() {
double pressure = 10; // bara
double temperature = 20; // C
double massFlowRate = 0.25; // kg/s
double heatTransferCoeff = 755; // W/m2K
double constantSurfaceTemperature = 100; // C

neqsim.thermo.system.SystemInterface testSystem = new neqsim.thermo.system.SystemSrkEos(
(273.15 + 45), ThermodynamicConstantsInterface.referencePressure);

testSystem.addComponent("water", 1.0);
testSystem.setMixingRule(2);
testSystem.init(0);
testSystem.useVolumeCorrection(true);
testSystem.setPressure(pressure, "bara");
testSystem.setTemperature(temperature, "C");
testSystem.setTotalFlowRate(massFlowRate, "kg/sec");

ThermodynamicOperations testOps = new ThermodynamicOperations(testSystem);
testOps.TPflash();
testSystem.initPhysicalProperties();

Stream stream_1 = new Stream("Stream1", testSystem);
stream_1.setFlowRate(massFlowRate, "kg/sec");

PipeBeggsAndBrills pipe = new PipeBeggsAndBrills("beggs and brils pipe 1", stream_1);
pipe.setPipeWallRoughness(0);
pipe.setLength(6.0);
pipe.setAngle(0);
pipe.setDiameter(0.05);
pipe.setNumberOfIncrements(1);
pipe.setConstantSurfaceTemperature(constantSurfaceTemperature, "C");
pipe.setHeatTransferCoefficient(heatTransferCoeff);
pipe.setNumberOfIncrements(10);

PipeBeggsAndBrills pipe2 = new PipeBeggsAndBrills("beggs and brils pipe 2", stream_1);
pipe2.setPipeWallRoughness(0);
pipe2.setLength(6.0);
pipe2.setAngle(0);
pipe2.setDiameter(0.05);
pipe2.setNumberOfIncrements(1);
pipe2.setConstantSurfaceTemperature(constantSurfaceTemperature, "C");


neqsim.process.processmodel.ProcessSystem operations =
new neqsim.process.processmodel.ProcessSystem();
operations.add(stream_1);
operations.add(pipe);
operations.add(pipe2);
operations.run();

double pressureOut = pipe.getOutletPressure();
double temperatureOut = pipe.getOutletTemperature() - 273.15;
double temperatureOut2 = pipe2.getOutletTemperature() - 273.15;
Assertions.assertEquals(temperatureOut, 57, 5);
Assertions.assertEquals(temperatureOut2, 40, 5);
}

}

0 comments on commit ec372df

Please sign in to comment.