Skip to content

Commit b11ee04

Browse files
DanWBRDanWBR
committed
feat: the bubble-point column solver relaxes its temperature step to converge harder columns
- when a bubble-point solve does not converge, retry with a smaller temperature step (0.25, then 0.1) on a reduced iteration budget before the wide-boiling fallback, so a cold-feed stabiliser and similar columns converge from an ordinary starting estimate (issue #73) - add Column.TemperatureStepFraction (default 0.5, the historical half step), used by both Wang-Henke solvers, on the Avalonia column editors and the fluent builders, so a user can preset the step and skip the failed attempts - add a cold naphtha stabiliser regression test and a setting round-trip test, and note the change in the release notes
1 parent 5290bb9 commit b11ee04

12 files changed

Lines changed: 168 additions & 55 deletions

File tree

‎content/readme.txt‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -44,6 +44,7 @@ Version 10.2.7
4444
- [CHG] The MCP server accepts a --python-path option, and a DWSIM_PYTHON_PATH environment variable for containers, so a Python Script unit operation can run on a headless server. The Python-path error message now names the distribution folder on Windows and the shared library file on Linux and macOS.
4545
- [CHG] The /api/check endpoint reports the open simulation's name and file path.
4646
- [FIX] A rigorous distillation or absorption column using the bubble-point method no longer reports a wrong converged result on a wide-boiling mixture: it now checks that the stage compositions have settled, not only the temperature and the vapour flow.
47+
- [CHG] A rigorous column solved by the bubble-point method now retries with a smaller temperature step when it does not converge, which solves cold-feed stabilisers and similar columns that failed from an ordinary starting estimate. The step is also a column setting, Temperature Step Fraction (default 0.5), on the editor and the fluent builders.
4748
- [FIX] A Recycle set to Global Convergence (Broyden) no longer corrupts its tear stream on the second solver pass, where the temperature, pressure and flow collapsed and the flow could go negative.
4849
- [FIX] A vertical gas-liquid separator's height now accounts for the liquid residence time, so changing the residence time changes the calculated size.
4950
- [FIX] Reopening a simulation saved by the cross-platform interface no longer shows blank panels (issue #72).

‎content/whatsnew.txt‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -3,6 +3,7 @@
33
- [CHG] The MCP server accepts a --python-path option, and a DWSIM_PYTHON_PATH environment variable for containers, so a Python Script unit operation can run on a headless server. The Python-path error message now names the distribution folder on Windows and the shared library file on Linux and macOS.
44
- [CHG] The /api/check endpoint reports the open simulation's name and file path.
55
- [FIX] A rigorous distillation or absorption column using the bubble-point method no longer reports a wrong converged result on a wide-boiling mixture: it now checks that the stage compositions have settled, not only the temperature and the vapour flow.
6+
- [CHG] A rigorous column solved by the bubble-point method now retries with a smaller temperature step when it does not converge, which solves cold-feed stabilisers and similar columns that failed from an ordinary starting estimate. The step is also a column setting, Temperature Step Fraction (default 0.5), on the editor and the fluent builders.
67
- [FIX] A Recycle set to Global Convergence (Broyden) no longer corrupts its tear stream on the second solver pass, where the temperature, pressure and flow collapsed and the flow could go negative.
78
- [FIX] A vertical gas-liquid separator's height now accounts for the liquid residence time, so changing the residence time changes the calculated size.
89
- [FIX] Reopening a simulation saved by the cross-platform interface no longer shows blank panels (issue #72).

‎engine/DWSIM.FluentAPI/Builders/AbsorptionColumnBuilder.cs‎

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -13,6 +13,8 @@ internal AbsorptionColumnBuilder(Flowsheet f, AbsorptionColumn o) : base(f, o) {
1313
public AbsorptionColumnBuilder WithTopPressure(Quantity p) { Object.SetTopPressure(p.SI); return this; }
1414
/// <summary>Sets <c>Column Pressure Drop</c> (SI) and returns this builder for chaining.</summary>
1515
public AbsorptionColumnBuilder WithColumnPressureDrop(Quantity dp) { Object.ColumnPressureDrop = dp.SI; return this; }
16+
/// <summary>Sets <c>Temperature Step Fraction</c> of the bubble-point solvers (0 to 1, default 0.5) and returns this builder for chaining.</summary>
17+
public AbsorptionColumnBuilder WithTemperatureStepFraction(double fraction) { Object.TemperatureStepFraction = fraction; return this; }
1618

1719
/// <summary>Sets <c>Feed</c> and returns this builder for chaining.</summary>
1820
public AbsorptionColumnBuilder WithFeed(MaterialStreamBuilder feed, int stageNumber)

‎engine/DWSIM.FluentAPI/Builders/DistillationColumnBuilder.cs‎

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -14,6 +14,8 @@ internal DistillationColumnBuilder(Flowsheet f, DistillationColumn o) : base(f,
1414
public DistillationColumnBuilder WithTopPressure(Quantity p) { Object.SetTopPressure(p.SI); return this; }
1515
/// <summary>Sets <c>Column Pressure Drop</c> (SI) and returns this builder for chaining.</summary>
1616
public DistillationColumnBuilder WithColumnPressureDrop(Quantity dp) { Object.ColumnPressureDrop = dp.SI; return this; }
17+
/// <summary>Sets <c>Temperature Step Fraction</c> of the bubble-point solvers (0 to 1, default 0.5) and returns this builder for chaining.</summary>
18+
public DistillationColumnBuilder WithTemperatureStepFraction(double fraction) { Object.TemperatureStepFraction = fraction; return this; }
1719

1820
/// <summary>Sets <c>Feed</c> and returns this builder for chaining.</summary>
1921
public DistillationColumnBuilder WithFeed(MaterialStreamBuilder feed, int stageNumber)

‎engine/DWSIM.FluentAPI/docs/api/unit-operations/columns.md‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -36,6 +36,7 @@ stage-by-stage rigorous solver.
3636
| `WithReboilerSpec(specType, value, units, compound = "")` | E.g. `"Product Molar Flow Rate"`. |
3737
| `WithTopPressure(p)` | Top-stage pressure. |
3838
| `WithColumnPressureDrop(dp)` | Total drop across the column. |
39+
| `WithTemperatureStepFraction(f)` | Fraction of the bubble-point stage-temperature update the Wang-Henke solvers apply each sweep, 0 to 1 (default 0.5). |
3940

4041
```csharp
4142
fs.AddDistillationColumn("T-101")

‎engine/DWSIM.UnitOperations/UnitOperations/RigorousColumn.vb‎

Lines changed: 9 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1893,6 +1893,15 @@ Namespace UnitOperations
18931893

18941894
Public Property ColumnPressureDrop As Double = Double.NaN
18951895

1896+
''' <summary>
1897+
''' Fraction of the bubble-point temperature update the Wang-Henke solvers apply on each sweep (0 to 1].
1898+
''' Each stage temperature moves from its previous value towards the newly computed bubble point by this fraction.
1899+
''' 0.5 is the historical behaviour and the right step for most columns. When the iteration does not converge
1900+
''' the Wang-Henke solver halves the step on its own, down to 0.1, before it tries the wide-boiling path, so a
1901+
''' column that needs a smaller step gets one without this being set; setting it skips those failed attempts.
1902+
''' </summary>
1903+
Public Property TemperatureStepFraction As Double = 0.5
1904+
18961905
Public Property TraySpacing As Double = 0.5 'm
18971906

18981907
Public Property EstimatedDiameter As Double = Double.NaN 'm

‎engine/DWSIM.UnitOperations/UnitOperations/RigorousColumnSolvers/BubblePoint.vb‎

Lines changed: 70 additions & 54 deletions
Original file line numberDiff line numberDiff line change
@@ -48,6 +48,10 @@ Namespace UnitOperations.Auxiliary.SepOps.SolvingMethods
4848
End Get
4949
End Property
5050

51+
' The temperature step in force for the current solve: the column's TemperatureStepFraction to
52+
' start with, lowered by SolveWithFallbacks when the iteration does not converge.
53+
Private _trelax As Double = 0.5
54+
5155
Public Function Solve(ByVal rc As Column, ByVal nc As Integer, ByVal ns As Integer, ByVal maxits As Integer,
5256
ByVal tol As Double(), ByVal F As Double(), ByVal V As Double(),
5357
ByVal Q As Double(), ByVal L As Double(),
@@ -129,7 +133,8 @@ Namespace UnitOperations.Auxiliary.SepOps.SolvingMethods
129133
Dim ObjFunctionValues As New List(Of Double)
130134
Dim ResultsVector As New List(Of Object)
131135

132-
Dim altmode As Boolean = False
136+
_trelax = 0.5
137+
If rc.TemperatureStepFraction > 0.0 AndAlso rc.TemperatureStepFraction <= 1.0 Then _trelax = rc.TemperatureStepFraction
133138

134139
If Not specC_OK And Not specR_OK Then
135140

@@ -354,20 +359,9 @@ Namespace UnitOperations.Auxiliary.SepOps.SolvingMethods
354359

355360
Dim result As Object = Nothing
356361

357-
Try
358-
result = Solve_Internal(rc, nc, ns, maxits, tolerance, F, V, Q, L, VSS, LSS, Kval,
359-
x, y, z, fc, HF, T, P, condt, stopatitnumber, eff,
360-
coltype, pp, newspecs, IdealK, IdealH, 0, flashalgs)
361-
altmode = False
362-
Catch ex As Exception
363-
altmode = True
364-
End Try
365-
366-
If altmode Then
367-
result = Solve_Internal(rc, nc, ns, maxits, tolerance, F, V, Q, L, VSS, LSS, Kval,
368-
x, y, z, fc, HF, T, P, condt, stopatitnumber, eff,
369-
coltype, pp, newspecs, IdealK, IdealH, 1, flashalgs)
370-
End If
362+
result = SolveWithFallbacks(rc, nc, ns, maxits, tolerance, F, V, Q, L, VSS, LSS, Kval,
363+
x, y, z, fc, HF, T, P, condt, stopatitnumber, eff,
364+
coltype, pp, newspecs, IdealK, IdealH, flashalgs)
371365

372366
T = result(0)
373367
V = result(1)
@@ -390,20 +384,9 @@ Namespace UnitOperations.Auxiliary.SepOps.SolvingMethods
390384

391385
If cspec.SpecValue < 0 Then cspec.SpecValue = -cspec.SpecValue
392386

393-
Try
394-
result = Solve_Internal(rc, nc, ns, maxits, tolerance, F, V, Q, L, VSS, LSS, Kval,
395-
x, y, z, fc, HF, T, P, condt, stopatitnumber, eff,
396-
coltype, pp, newspecs, IdealK, IdealH, 0, flashalgs)
397-
altmode = False
398-
Catch ex As Exception
399-
altmode = True
400-
End Try
401-
402-
If altmode Then
403-
result = Solve_Internal(rc, nc, ns, maxits, tolerance, F, V, Q, L, VSS, LSS, Kval,
404-
x, y, z, fc, HF, T, P, condt, stopatitnumber, eff,
405-
coltype, pp, newspecs, IdealK, IdealH, 1, flashalgs)
406-
End If
387+
result = SolveWithFallbacks(rc, nc, ns, maxits, tolerance, F, V, Q, L, VSS, LSS, Kval,
388+
x, y, z, fc, HF, T, P, condt, stopatitnumber, eff,
389+
coltype, pp, newspecs, IdealK, IdealH, flashalgs)
407390

408391
'Return New Object() {Tj, Vj, Lj, VSSj, LSSj, yc, xc, K, Q, ic, t_error}
409392

@@ -529,20 +512,9 @@ Namespace UnitOperations.Auxiliary.SepOps.SolvingMethods
529512

530513
Dim result As Object = Nothing
531514

532-
Try
533-
result = Solve_Internal(rc, nc, ns, maxits, tolerance, F, V, Q, L, VSS, LSS, Kval,
534-
x, y, z, fc, HF, T, P, condt, stopatitnumber, eff,
535-
coltype, pp, newspecs, IdealK, IdealH, 0, flashalgs)
536-
altmode = False
537-
Catch ex As Exception
538-
altmode = True
539-
End Try
540-
541-
If altmode Then
542-
result = Solve_Internal(rc, nc, ns, maxits, tolerance, F, V, Q, L, VSS, LSS, Kval,
543-
x, y, z, fc, HF, T, P, condt, stopatitnumber, eff,
544-
coltype, pp, newspecs, IdealK, IdealH, 1, flashalgs)
545-
End If
515+
result = SolveWithFallbacks(rc, nc, ns, maxits, tolerance, F, V, Q, L, VSS, LSS, Kval,
516+
x, y, z, fc, HF, T, P, condt, stopatitnumber, eff,
517+
coltype, pp, newspecs, IdealK, IdealH, flashalgs)
546518

547519
T = result(0)
548520
V = result(1)
@@ -743,20 +715,59 @@ Namespace UnitOperations.Auxiliary.SepOps.SolvingMethods
743715

744716
Else
745717

718+
Return SolveWithFallbacks(rc, nc, ns, maxits, tolerance, F, V, Q, L, VSS, LSS, Kval, x, y, z, fc, HF, T, P, condt,
719+
stopatitnumber, eff, coltype, pp, specs, IdealK, IdealH, flashalgs)
720+
721+
End If
722+
723+
End Function
724+
725+
''' <summary>
726+
''' One bubble-point solve with its fallbacks. The column's temperature step first, at the full
727+
''' iteration budget, so a column that converges today is unaffected. When that does not converge,
728+
''' the step is halved down to 0.1 (0.25, then 0.1) and the solve repeated from the same start, on
729+
''' a reduced budget so a column that no step converges does not pay the full budget several times
730+
''' over. When the fallbacks fail too, the wide-boiling path (Mode 1) as before. The step that
731+
''' worked is kept for the rest of this solve, since the specification loops re-solve many times.
732+
''' </summary>
733+
Private Function SolveWithFallbacks(rc As Column, nc As Integer, ns As Integer, maxits As Integer,
734+
tolerance As Double, F As Double(), V As Double(),
735+
Q As Double(), L As Double(),
736+
VSS As Double(), LSS As Double(), Kval()() As Double,
737+
x()() As Double, y()() As Double, z()() As Double,
738+
fc()() As Double,
739+
HF As Double(), T As Double(), P As Double(),
740+
condt As DistillationColumn.condtype,
741+
stopatitnumber As Integer,
742+
eff() As Double,
743+
coltype As Column.ColType,
744+
pp As PropertyPackages.PropertyPackage,
745+
specs As Dictionary(Of String, SepOps.ColumnSpec),
746+
IdealK As Boolean, IdealH As Boolean,
747+
flashalgs As List(Of FlashAlgorithm)) As Object
748+
749+
' First attempt at the column's own step (kept from a previous solve if it was lowered),
750+
' at the full budget: a column that converges at its step is not slowed down.
751+
Try
752+
Return Solve_Internal(rc, nc, ns, maxits, tolerance, F, V, Q, L, VSS, LSS, Kval, x, y, z, fc, HF, T, P, condt, stopatitnumber, eff, coltype, pp, specs, IdealK, IdealH, 0, flashalgs)
753+
Catch ex As Exception
754+
End Try
755+
756+
' Exploratory rungs at a smaller step, on a reduced budget so a hard column does not pay the
757+
' full budget on every rung. Enough to converge the columns these rungs are for.
758+
Dim fbmaxits As Integer = Math.Max(75, maxits \ 2)
759+
Do
760+
If _trelax <= 0.1 Then Exit Do
761+
_trelax = If(_trelax > 0.25, 0.25, 0.1)
762+
pp.Flowsheet?.ShowMessage(rc.GraphicObject.Tag + ": [BP Solver] did not converge, retrying with the temperature step relaxed to " + _trelax.ToString("0.00"), IFlowsheet.MessageType.Information)
746763
Try
747-
Return Solve_Internal(rc, nc, ns, maxits, tolerance, F, V, Q, L, VSS, LSS, Kval, x, y, z, fc, HF, T, P, condt,
748-
stopatitnumber, eff, coltype, pp, specs, IdealK, IdealH, 0, flashalgs)
749-
altmode = False
764+
Return Solve_Internal(rc, nc, ns, fbmaxits, tolerance, F, V, Q, L, VSS, LSS, Kval, x, y, z, fc, HF, T, P, condt, stopatitnumber, eff, coltype, pp, specs, IdealK, IdealH, 0, flashalgs)
750765
Catch ex As Exception
751-
altmode = True
752766
End Try
767+
Loop
753768

754-
If altmode Then
755-
Return Solve_Internal(rc, nc, ns, maxits, tolerance, F, V, Q, L, VSS, LSS, Kval, x, y, z, fc, HF, T, P, condt,
756-
stopatitnumber, eff, coltype, pp, specs, IdealK, IdealH, 1, flashalgs)
757-
End If
758-
759-
End If
769+
' Wide-boiling path as the last resort, at the full budget as before.
770+
Return Solve_Internal(rc, nc, ns, maxits, tolerance, F, V, Q, L, VSS, LSS, Kval, x, y, z, fc, HF, T, P, condt, stopatitnumber, eff, coltype, pp, specs, IdealK, IdealH, 1, flashalgs)
760771

761772
End Function
762773

@@ -858,6 +869,11 @@ Namespace UnitOperations.Auxiliary.SepOps.SolvingMethods
858869

859870
Dim ic As Integer
860871
Dim t_error, t_error_ant, vf_error, x_error, xcerror(ns) As Double
872+
' Under-relaxation of the bubble-point temperature update: the step in force for this
873+
' solve (the column's TemperatureStepFraction, lowered by SolveWithFallbacks on a failure).
874+
' 0.5 is the historical half step. Anything outside (0, 1] falls back to it.
875+
Dim trelax As Double = _trelax
876+
If trelax <= 0.0 OrElse trelax > 1.0 Then trelax = 0.5
861877
Dim Tj(ns), Tj_ant(ns), dTj(ns) As Double
862878
Dim Fj(ns), Lj(ns), Vj(ns), Vj_ant(ns), dVj(ns), xc(ns)(), xc0(ns)(), fcj(ns)(), yc(ns)(), lc(ns)(), vc(ns)(), zc(ns)(), K(ns)(), Kant(ns)() As Double
863879
Dim Hfj(ns), Hv(ns), Hl(ns) As Double
@@ -1327,7 +1343,7 @@ Namespace UnitOperations.Auxiliary.SepOps.SolvingMethods
13271343
Else
13281344
'If ic < 5 Or ic > 100 Then
13291345
For i = 0 To ns
1330-
Tj(i) = Tj(i) / 2 + Tj_ant(i) / 2
1346+
Tj(i) = trelax * Tj(i) + (1.0 - trelax) * Tj_ant(i)
13311347
Next
13321348
'End If
13331349
End If

‎engine/DWSIM.UnitOperations/UnitOperations/RigorousColumnSolvers/BubblePoint2.vb‎

Lines changed: 7 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1314,6 +1314,12 @@ Namespace UnitOperations.Auxiliary.SepOps.SolvingMethods
13141314
Dim doparallel As Boolean = Settings.EnableParallelProcessing
13151315

13161316
Dim ic As Integer
1317+
' Under-relaxation of the bubble-point temperature update, a column setting
1318+
' (default 0.5, the historical half step). Anything outside (0, 1] falls back to it.
1319+
Dim trelax As Double = 0.5
1320+
If rc.TemperatureStepFraction > 0.0 AndAlso rc.TemperatureStepFraction <= 1.0 Then
1321+
trelax = rc.TemperatureStepFraction
1322+
End If
13171323
Dim t_error, t_error_ant, vf_error, xcerror(ns) As Double
13181324
Dim Tj(ns), Tj_ant(ns), dTj(ns) As Double
13191325
Dim Fj(ns), Lj(ns), Vj(ns), Vj_ant(ns), dVj(ns), xc(ns)(), xc0(ns)(), fcj(ns)(), yc(ns)(), lc(ns)(), vc(ns)(), zc(ns)(), K(ns)(), Kant(ns)() As Double
@@ -1765,7 +1771,7 @@ Namespace UnitOperations.Auxiliary.SepOps.SolvingMethods
17651771
Next
17661772
Else
17671773
For i = 0 To ns
1768-
Tj(i) = Tj(i) / 2 + Tj_ant(i) / 2
1774+
Tj(i) = trelax * Tj(i) + (1.0 - trelax) * Tj_ant(i)
17691775
Next
17701776
End If
17711777

0 commit comments

Comments
 (0)