--- title: "Original Source Code" author: "Bill Denney" date: "`r Sys.Date()`" output: rmarkdown::html_vignette vignette: > %\VignetteIndexEntry{Original Source Code} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- ```{r setup, include = FALSE} knitr::opts_chunk$set( collapse = TRUE, comment = "#>" ) ``` Below are the original Visual Basic source code used for the translation to R. It is included to allow users to inspect the fidelity of the code translation. This program is intended for research use, only. The code within is translated from Visual Basic code based on Werness, et al 1985 to R. The Visual Basic code was kindly provided by Dr. John Lieske of the Mayo Clinic. # Reference Werness PG, Brown CM, Smith LH, Finlayson B. Equil2: A Basic Computer Program for the Calculation of Urinary Saturation. Journal of Urology. 1985;134(6):1242-1244. doi:10.1016/S0022-5347(17)47703-2 # EQUIL-2 The code below is the basis of the `equil2()` function. ``` VERSION 5.00 Begin VB.Form Form1 Caption = "Form1" ClientHeight = 8355 ClientLeft = 60 ClientTop = 345 ClientWidth = 5100 LinkTopic = "Form1" ScaleHeight = 8355 ScaleWidth = 5100 StartUpPosition = 3 'Windows Default Begin VB.TextBox Text15 Height = 375 Left = 1680 TabIndex = 2 Top = 480 Width = 3375 End Begin VB.TextBox Text14 Height = 375 Left = 1680 TabIndex = 1 Top = 0 Width = 3375 End Begin VB.TextBox Text11 Height = 375 Left = 1680 TabIndex = 13 Top = 6840 Width = 1215 End Begin VB.TextBox Text10 Height = 375 Left = 1680 TabIndex = 12 Top = 6360 Width = 1215 End Begin VB.TextBox Text9 Height = 375 Left = 1680 TabIndex = 11 Top = 5880 Width = 1215 End Begin VB.TextBox Text8 Height = 375 Left = 1680 TabIndex = 10 Top = 5400 Width = 1215 End Begin VB.TextBox Text7 Height = 375 Left = 1680 TabIndex = 9 Top = 4920 Width = 1215 End Begin VB.TextBox Text6 Height = 375 Left = 1680 TabIndex = 8 Top = 4440 Width = 1215 End Begin VB.TextBox Text5 Height = 375 Left = 1680 TabIndex = 7 Top = 3960 Width = 1215 End Begin VB.TextBox Text4 Height = 375 Left = 1680 TabIndex = 6 Top = 3480 Width = 1215 End Begin VB.TextBox Text3 Height = 375 Left = 1680 TabIndex = 5 Top = 3000 Width = 1215 End Begin VB.TextBox Text2 Height = 375 Left = 1680 TabIndex = 4 Top = 2520 Width = 1215 End Begin VB.TextBox Text1 Height = 375 Left = 1680 TabIndex = 3 Top = 2040 Width = 1215 End Begin VB.CommandButton Command1 Caption = "Calculate now" Height = 615 Left = 3240 TabIndex = 0 Top = 7800 Width = 1935 End Begin VB.Label Label1 Caption = "Enter Collection Date, Time" Height = 375 Index = 14 Left = 120 TabIndex = 26 Top = 480 Width = 1455 End Begin VB.Label Label1 Caption = "Enter Patient Name, Clinic Number" Height = 495 Index = 13 Left = 120 TabIndex = 25 Top = 0 Width = 1455 End Begin VB.Label Label1 Caption = "pH" Height = 255 Index = 11 Left = 0 TabIndex = 24 Top = 6480 Width = 1455 End Begin VB.Label Label1 Caption = "Uric Acid (mg/dl)" Height = 375 Index = 10 Left = 0 TabIndex = 23 Top = 6840 Width = 1455 End Begin VB.Label Label1 Caption = "Citrate (mg/dl)" Height = 375 Index = 9 Left = 0 TabIndex = 22 Top = 5520 Width = 1455 End Begin VB.Label Label1 Caption = "Oxalate (mg/dl)" Height = 375 Index = 8 Left = 0 TabIndex = 21 Top = 6000 Width = 1455 End Begin VB.Label Label1 Caption = "Phosphate (mg/dl)" Height = 375 Index = 7 Left = 0 TabIndex = 20 Top = 4560 Width = 1455 End Begin VB.Label Label1 Caption = "Sulfate (mg/dl)" Height = 375 Index = 6 Left = 0 TabIndex = 19 Top = 5040 Width = 1455 End Begin VB.Label Label1 Caption = "Chloride (Meq/L)" Height = 375 Index = 5 Left = 0 TabIndex = 18 Top = 4080 Width = 1455 End Begin VB.Label Label1 Caption = "Sodium (Meq/L)" Height = 375 Index = 3 Left = 0 TabIndex = 17 Top = 2160 Width = 1455 End Begin VB.Label Label1 Caption = "Potassium (Meq/L)" Height = 375 Index = 2 Left = 0 TabIndex = 16 Top = 2640 Width = 1455 End Begin VB.Label Label1 Caption = "Magnesium (mg/dl)" Height = 375 Index = 1 Left = 0 TabIndex = 15 Top = 3600 Width = 1455 End Begin VB.Label Label1 Caption = "Calcium (mg/dl)" Height = 375 Index = 0 Left = 0 TabIndex = 14 Top = 3120 Width = 1455 End End Attribute VB_Name = "Form1" Attribute VB_GlobalNameSpace = False Attribute VB_Creatable = False Attribute VB_PredeclaredId = True Attribute VB_Exposed = False Dim p As String Dim q As String Dim a(120) As Double Sub Command1_Click() Call calc End Sub Sub Form_Load() Debug.Print "hello" End Sub Sub calc() Rem Open "ssinput.txt" For Input As #1 Debug.Print "Now in calc routine." 10 Dim S(65) 20 S(1) = 1730000000000# 30 S(2) = 14900000# 40 S(3) = 162 50 S(4) = 145.5 60 S(5) = 20750 70 S(6) = 2640000! 80 S(7) = 55210! 90 S(8) = 1247 100 S(9) = 278000! 110 S(10) = 12.9 120 S(11) = 5.433 130 S(12) = 13.4 140 S(13) = 8.5 150 S(14) = 216 160 S(15) = 33.1 170 S(16) = 251.2 180 S(17) = 10 190 S(18) = 8.831001 200 S(19) = 13.4 210 S(20) = 12.6 220 S(21) = 143 230 S(22) = 3597000! 240 S(23) = 685 250 S(24) = 31.3 260 S(25) = 229.6 270 S(26) = 2746 280 S(27) = 17.3 290 S(28) = 71.4 300 S(29) = 60000! 310 S(30) = 505.2 320 S(31) = 12.5 330 S(32) = 3460000! 340 S(33) = 1014 350 S(34) = 31.9 360 S(35) = 188.4 370 S(36) = 4020 380 S(37) = 4.75 390 S(38) = 5.93 400 S(39) = 69900! 410 S(40) = 316.7 420 S(41) = 5 430 S(42) = 10 440 S(43) = 12.9 450 S(44) = 13 460 S(45) = 8.5 470 S(46) = 58800000# 480 S(47) = 2440000000# 490 S(48) = 4970000! 500 S(49) = 171 510 S(50) = 7.05 520 S(51) = 562000! 530 S(52) = 5500 540 S(53) = 794000000# 550 S(54) = 23.1 560 S(55) = 380.19 570 S(56) = 19770000# 580 S(57) = 1995000000# 590 S(58) = 55.6 600 S(59) = 1940000! 610 S(60) = 16600000000# 620 S(61) = 0.00123 630 S(62) = 18.6 640 S(63) = 1.03 650 S(64) = 1897 660 S(65) = 946 670 Rem Dim a(120) 680 Dim T(11) 690 Rem Open "LPT1:" For Output As #1 691 Rem CLS: Key OFF: Color 7, 1, 0 695 Rem Input #1, z$ 697 Rem Input #1, z$ 700 Rem Input #1, a(120) 701 a(1) = a(1) / 1000 702 Rem Print "K";: Input #1, a(120) 703 a(2) = a(2) / 1000 704 Rem Print "CA";: Input #1, a(120) 705 a(3) = a(3) / 4008 706 Rem Print "MG";: Input #1, a(120) 707 a(4) = a(4) / 2431 708 Rem Print "PYRO";: Input #1, a(11) 709 Rem a(11)=a(11)/2431 710 Rem Print "NH4";: Input #1, a(120) 711 a(5) = a(5) / 1000 712 Rem Print "CL";: Input #1, a(120) 713 a(91) = a(91) / 1000 714 Rem Print "CO2";: Input #1, a(120) 715 Rem a(31)=a(31)/4401 716 Rem Print "P";: Input #1, a(120) 717 a(6) = a(6) / 3097 718 Rem Print "S";: Input #1, a(120) 719 a(7) = a(7) / 3206 720 Rem Print "CIT";: Input #1, a(120) 721 a(9) = a(9) / 19212 722 Rem Print "OX";: Input #1, a(120) 723 a(8) = a(8) / 8802 724 Rem Print "PH";: Input #1, a(120) 726 Rem Print "UR";: Input #1, a(120) 727 a(10) = a(10) / 16800 800 Rem Close #1 820 a(26) = 10 ^ (-a(24)) 850 F1 = 0.7 860 F2 = 0.3 870 F3 = 0.1 880 F4 = 0.02 890 O0 = 0 900 O1 = 0 910 O2 = 0 920 O3 = 0 930 O4 = 0 940 For I = 1 To 10 950 a(12 + I) = 0.1 * a(I) 960 Next I 970 a(12) = 0.1 * a(11) 980 a(23) = 0.01 * a(12) 990 For I = 1 To 50 1000 a(25) = 10 ^ (-13.593 + a(24)) 1010 a(27) = S(1) * a(26) * a(18) * F3 / F2 1020 a(28) = S(2) * a(26) * a(27) * F2 / F1 1030 a(29) = S(3) * a(26) * a(28) * F1 1040 a(30) = S(61) * a(31) * S(58) 1050 a(32) = a(30) / (a(26) * S(59) * F1) 1060 a(33) = a(32) * F1 / (a(26) * S(60) * F2) 1070 a(34) = S(62) * a(13) * a(33) * F2 1080 a(35) = S(63) * a(13) * a(34) * F1 * F1 1090 a(36) = S(64) * a(15) * a(33) * F2 * F2 1100 a(37) = S(65) * a(16) * a(33) * F2 * F2 1110 a(38) = S(4) * a(26) * a(19) * F2 / F1 1120 a(39) = S(5) * a(26) * a(20) * F2 / F1 1130 a(40) = S(6) * a(26) * a(21) * F3 / F2 1140 a(41) = S(7) * a(26) * a(40) * F2 / F1 1150 a(42) = S(8) * a(26) * a(41) * F1 1160 a(43) = S(9) * a(26) * a(22) * F1 1170 a(44) = S(10) * a(13) * a(27) * F2 1180 a(45) = S(11) * a(13) * a(19) * F2 1190 a(46) = S(12) * a(13) * a(20) * F2 1200 a(47) = S(13) * a(13) * a(21) * F3 * F1 / F2 1210 a(48) = S(14) * a(13) * a(12) * F1 * F4 / F3 1220 a(49) = S(16) * a(13) * a(48) * F1 * F3 / F2 1230 a(50) = S(15) * a(13) * a(23) * F1 * F3 / F2 1240 a(51) = S(17) * a(14) * a(27) * F2 1250 a(52) = S(18) * a(14) * a(19) * F2 1260 a(53) = S(19) * a(14) * a(20) * F2 1270 a(54) = S(20) * a(14) * a(21) * F3 * F1 / F2 1280 a(55) = S(21) * a(14) * a(12) * F1 * F4 / F3 1290 a(56) = S(22) * a(15) * a(18) * F3 * F2 / F1 1300 a(57) = S(23) * a(15) * a(27) * F2 * F2 1310 a(58) = S(24) * a(15) * a(28) * F2 1320 a(59) = S(25) * a(15) * a(19) * F2 * F2 1330 a(60) = S(26) * a(15) * a(20) * F2 * F2 1340 a(61) = S(28) * a(15) * a(60) 1350 a(62) = S(27) * a(60) * a(20) 1360 a(63) = S(29) * a(15) * a(21) * F3 * F2 / F1 1370 a(64) = S(30) * a(15) * a(40) * F2 * F2 1380 a(65) = S(31) * a(15) * a(41) * F2 1390 a(66) = S(32) * a(16) * a(18) * F3 * F2 / F1 1400 a(67) = S(33) * a(16) * a(27) * F2 * F2 1410 a(68) = S(34) * a(16) * a(28) * F2 1420 a(69) = S(35) * a(16) * a(19) * F2 * F2 1430 a(70) = S(36) * a(16) * a(20) * F2 * F2 1440 a(71) = S(37) * a(16) * a(70) 1450 a(72) = S(38) * a(70) * a(20) 1460 a(73) = S(39) * a(16) * a(21) * F3 * F2 / F1 1470 a(74) = S(40) * a(16) * a(40) * F2 * F2 1480 a(75) = S(41) * a(16) * a(41) * F2 1490 a(76) = S(42) * a(17) * a(27) * F2 1500 a(77) = S(43) * a(17) * a(19) * F2 1510 a(78) = S(44) * a(17) * a(20) * F2 1520 a(79) = S(45) * a(17) * a(21) * F3 * F1 / F2 1540 a(23) = S(47) * a(26) * a(12) * F4 / F3 1550 a(81) = S(48) * a(26) * a(23) * F3 / F2 1560 a(82) = S(49) * a(26) * a(81) * F2 / F1 1570 a(83) = S(50) * a(26) * a(82) * F1 1580 a(84) = S(51) * a(15) * a(12) * F4 1590 a(85) = S(52) * a(15) * a(23) * F2 * F3 / F1 1600 a(86) = S(53) * a(15) * a(25) * a(12) * F4 * F2 / F3 1610 a(87) = S(54) * a(15) * a(25) * F2 / F1 1620 a(88) = S(55) * a(16) * a(25) * F2 / F1 1630 a(89) = S(56) * a(16) * a(12) * F4 1640 a(90) = S(57) * a(16) * a(25) * a(12) * F2 * F4 / F3 1650 T(0) = a(13) + a(44) + a(45) + a(46) + a(47) + a(48) + 2 * a(49) + a(50) + a(34) + 2 * a(35) 1660 T(1) = a(14) + a(51) + a(52) + a(53) + a(54) + a(55) 1670 T(2) = a(17) + a(76) + a(77) + a(78) + a(79) 1680 T(3) = a(15) + a(56) + a(57) + a(58) + a(59) + a(60) + 2 * a(61) + a(63) + a(64) + a(62) + a(65) + a(36) + a(85) + a(84) + a(86) + a(87) 1690 T(4) = a(16) + a(66) + a(67) + a(68) + a(69) + a(70) + 2 * a(71) + a(73) + a(74) + a(75) + a(72) + a(37) + a(88) + a(89) + a(90) 1700 T(5) = a(18) + a(27) + a(28) + a(29) + a(44) + a(51) + a(56) + a(57) + a(58) + a(66) + a(67) + a(68) + a(76) 1710 T(6) = a(19) + a(38) + a(45) + a(52) + a(59) + a(69) + a(77) 1720 T(7) = a(20) + a(39) + a(60) + a(61) + a(70) + a(71) + a(46) + a(53) + a(78) + 2 * a(62) + 2 * a(72) 1730 T(8) = a(21) + a(40) + a(41) + a(42) + a(47) + a(54) + a(79) + a(63) + a(64) + a(65) + a(73) + a(74) + a(75) 1740 T(9) = a(12) + a(23) + a(81) + a(82) + a(83) + a(84) + a(85) + a(86) + a(89) + a(90) + a(48) + a(49) + a(50) + a(55) 1750 T(10) = a(22) + a(43) 1760 T(11) = a(33) + a(32) + a(30) + a(34) + a(35) + a(36) + a(37) 1770 For I1 = O To 11 1780 If T(I1) = 0 Then T(I1) = 1E-20 1790 Next I1 1800 a(13) = a(1) * a(13) / T(0) 1810 a(14) = a(2) * a(14) / T(1) 1820 a(15) = a(3) * a(15) / T(3) 1830 a(16) = a(4) * a(16) / T(4) 1840 a(17) = a(5) * a(17) / T(2) 1850 a(18) = a(6) * a(18) / T(5) 1860 a(19) = a(7) * a(19) / T(6) 1870 a(20) = a(8) * a(20) / T(7) 1880 a(21) = a(9) * a(21) / T(8) 1890 a(22) = a(10) * a(22) / T(10) 1900 a(12) = a(11) * a(12) / T(9) 1910 a(33) = a(31) * a(33) / T(11) 1920 S1 = (a(25) + a(26)) / F1 + a(13) + a(14) + a(17) + a(22) + a(91) + a(44) + a(45) + a(46) + a(34) + a(51) + a(52) + a(53) + a(76) + a(77) 1930 S1 = S1 + a(78) + a(56) + a(58) + a(63) + a(65) + a(85) + a(87) + a(28) + a(32) + a(38) + a(39) + a(41) + a(82) + a(66) + a(68) + a(73) + a(75) + a(88) 1940 S2 = 4 * (a(15) + a(16) + a(19) + a(20) + a(33) + a(47) + a(49) + a(54) + a(79) + a(84) + a(50) + a(89) + a(27) + a(40) + a(81) + a(61) + a(71) + a(62) + a(72)) 1950 S3 = 9 * (a(18) + a(21) + a(48) + a(55) + a(23) + a(86) + a(90)) 1960 S4 = 16 * a(12) 1970 S5 = (S1 + S2 + S3 + S4) / 2 1980 If S5 > 1 Then S5 = 1 1990 If S5 < 0.000001 Then S5 = 0.000001 2000 S6 = Sqr(S5) 2010 F1 = Exp(-1.20218 * ((S6 / (1 + S6)) - 0.285 * S5)) 2020 F2 = F1 ^ 4 2030 F3 = F1 ^ 9 2040 F4 = F1 ^ 16 2050 If a(15) = 0 Then GoTo 2070 2060 If Abs((a(15) - O0) / a(15)) > 0.0001 Then GoTo 2160 2070 If a(16) = 0 Then GoTo 2090 2080 If Abs((a(16) - O1) / a(16)) > 0.0001 Then GoTo 2160 2090 If a(18) = 0 Then GoTo 2110 2100 If Abs((a(18) - O2) / a(18)) > 0.0001 Then GoTo 2160 2110 If a(20) = 0 Then GoTo 2130 2120 If Abs((a(20) - O3) / a(20)) > 0.0001 Then GoTo 2160 2130 If a(21) = 0 Then GoTo 2150 2140 If Abs((a(21) - O4) / a(21)) > 0.0001 Then GoTo 2160 2150 GoTo 2230 2160 O0 = a(15) 2170 O1 = a(16) 2180 O2 = a(18) 2190 O3 = a(20) 2200 O4 = a(21) 2210 Next I 2220 Print "50 CYCLES WITHOUT CONVERGE" 2230 a(92) = I 2231 a(93) = S5 2232 a(94) = F1 2234 a(95) = F2 2236 a(96) = F3 2238 a(97) = F4 2250 Open "ssout.txt" For Output As #1 2257 Print #1, "Patient: "; Tab(30); p$ 2258 Print #1, "Analysis Date :"; Tab(30); Date$ 2259 Print #1, "Collection Date and Time:"; Tab(30); q$: Print #1, 2260 Print #1, "Sodium", a(1) * 1000 2261 Rem Tab(43); "[NAHPP]"; Tab(57); A(50) 2270 Print #1, "Potassium", a(2) * 1000 2271 Rem Tab(43); "[KHPO4]"; Tab(57); A(51) 2280 Print #1, "Calcium", a(3) * 4008 2281 Rem Tab(43); "[KSO4]"; Tab(57); A(52) 2290 Print #1, "Magnesium", a(4) * 2431 2291 Rem Tab(43); "[KOX]"; Tab(57); A(53) 2300 Print #1, "Ammonia", a(5) * 1000 2301 Rem Tab(43); "[KCIT]"; Tab(57); A(54) 2302 Print #1, "Chloride", a(91) * 1000 2310 Print #1, "Phosphate", a(6) * 3097 2311 Rem Tab(43); "[KPP]"; Tab(57); A(55) 2320 Print #1, "Sulfate", a(7) * 3206 2321 Rem Tab(43); "[CAPO4]"; Tab(57); A(56) 2330 Print #1, "Oxalate", a(8) * 8802 2331 Rem Tab(43); "[CAHPO4]"; Tab(57); A(57) 2340 Print #1, "Citrate", a(9) * 19212 2341 Rem Tab(43); "[CAH2P04]"; Tab(57); A(58) 2342 Print #1, "pH", a(24) 2350 Print #1, "Urate", a(10) * 16800 2351 Rem Tab(43); "[CASO4]"; Tab(57); A(59) 2360 Rem Print #1, "PP", A(11) 2361 Rem Tab(43); "[CAOX]"; Tab(57); A(60) 2370 Rem Print #1, "[PP]", A(12), Tab(43); "[CA2OX]"; Tab(57); A(61) 2380 Rem Print #1, "[NA]", A(13), Tab(43); "[CAOX2]"; Tab(57); A(62) 2390 Rem Print #1, "[K]", A(14), Tab(43); "[CACIT]"; Tab(57); A(63) 2400 Rem Print #1, "[CA]", A(15), Tab(43); "[CAHCIT]"; Tab(57); A(64) 2410 Rem Print #1, "[MG]", A(16), Tab(43); "[CAH2CIT]"; Tab(57); A(65) 2420 Rem Print #1, "[NH4]", A(17), Tab(43); "[MGPO4]"; Tab(57); A(66) 2430 Rem Print #1, "[PO4]", A(18), Tab(43); "[MGHPO4]"; Tab(57); A(67) 2440 Rem Print #1, "[SO4]", A(19), Tab(43); "[MGH2PO4]"; Tab(57); A(68) 2450 Rem Print #1, "[OX]", A(20), Tab(43); "[MGSO4]"; Tab(57); A(69) 2460 Rem Print #1, "[CIT]", A(21), Tab(43); "[MGOX]"; Tab(57); A(70) 2470 Rem Print #1, "[HU]", A(22), Tab(43); "[MG2OX]"; Tab(57); A(71) 2480 Rem Print #1, "[HPP]", A(23), Tab(43); "[MGOX2]", Tab(57); A(72) 2490 Rem Print #1, "PH", A(24), Tab(43); "[MGCIT]"; Tab(57); A(73) 2500 Rem Print #1, "(OH)", A(25), Tab(43); "[MGHCIT]"; Tab(57); A(74) 2510 Rem Print #1, "(H)", A(26), Tab(43); "[MGH2CIT]"; Tab(57); A(75) 2520 Rem Print #1, "[HPO4]", A(27), Tab(43); "[NH4HPO4]"; Tab(57); A(76) 2530 Rem Print #1, "[H2PO4]", A(28), Tab(43); "[NH4SO4]"; Tab(57); A(77) 2540 Rem Print #1, "[H3PO4]", A(29), Tab(43); "[NH4OX]"; Tab(57); A(78) 2550 Rem Print #1, "[H2CO3]", A(30), Tab(43); "[NH4CIT]"; Tab(57); A(79) 2560 Rem Print #1, "CO2", A(31) 2570 Rem Print #1, "[HCO3]", A(32), Tab(43); "[H2PP]"; Tab(57); A(81) 2580 Rem Print #1, "[CO3]", A(33), Tab(43); "[H3PP]"; Tab(57); A(82) 2590 Rem Print #1, "[NACO3]", A(34), Tab(43); "[H4PP]"; Tab(57); A(83) 2600 Rem Print #1, "[NA2CO3]", A(35), Tab(43); "[CAPP]"; Tab(57); A(84) 2610 Rem Print #1, "[CACO3]", A(36), Tab(43); "[CAHPP]"; Tab(57); A(85) 2620 Rem Print #1, "[MGCO3]", A(37), Tab(43); "[CAOHPP]"; Tab(57); A(86) 2630 Rem Print #1, "[HSO4]", A(38), Tab(43); "[CAOH]"; Tab(57); A(87) 2640 Rem Print #1, "[HOX]", A(39), Tab(43); "[MGOH]"; Tab(57); A(88) 2650 Rem Print #1, "[HCIT]", A(40), Tab(43); "[MGPP]"; Tab(57); A(89) 2660 Rem Print #1, "[H2CIT]", A(41), Tab(43); "[MGOHPP]"; Tab(57); A(90) 2670 Rem Print #1, "[H3CIT]", A(42), Tab(43); "[CL]"; Tab(57); A(91) 2680 Rem Print #1, "[H2U]", A(43), Tab(43); "CYCLES"; Tab(57); A(92) 2690 Rem Print #1, "[NAHPO4]", A(44), Tab(43); "I.S."; Tab(57); A(93) 2700 Rem Print #1, "[NASO4]", A(45), Tab(43); "F1"; Tab(57); A(94) 2710 Rem Print #1, "[NAOX]", A(46), Tab(43); "F2"; Tab(57); A(95) 2720 Rem Print #1, "[NACIT]", A(47), Tab(43); "F3"; Tab(57); A(96) 2730 Rem Print #1, "[NAPP]", A(48), Tab(43); "F4"; Tab(57); A(97) 2740 Rem Print #1, "[NA2PP]", A(49) 2750 Print #1, "Cycles", a(92) 2760 Print #1, 2770 Print #1, Tab(20); "SS", Tab(43); "DG" 2780 a(100) = a(60) / 0.00000616 2790 a(101) = a(15) * a(27) * F2 * F2 / 0.000000237 2800 X1 = F2 * a(15) * 1000 2810 X2 = F3 * a(18) * 10000000000# 2820 a(102) = (X1 ^ 5) * (X2 ^ 3) * a(25) / 1.45E-14 2830 a(103) = F1 * F2 * F3 * a(16) * a(17) * a(18) / 0.000000000000115 2840 a(104) = a(43) / 0.000261 2850 a(105) = F1 * F1 * a(22) * a(13) / 0.0000279 2860 a(106) = F1 * F1 * a(22) * a(17) / 0.000036 2870 a(107) = F1 * F1 * a(22) * a(14) / 0.0000963 2875 If a(100) <= 0 Then a(108) = 0: GoTo 2885 2880 a(108) = 1.2935 * Log(a(100)) 2885 If a(101) <= 0 Then a(109) = 0: GoTo 2895 2890 a(109) = 1.2935 * Log(a(101)) 2895 If a(102) <= 0 Then a(110) = 0: GoTo 2905 2900 a(110) = 0.28744 * Log(a(102)) 2905 If a(103) <= 0 Then a(111) = 0: GoTo 2915 2910 a(111) = 0.8623 * Log(a(103)) 2915 If a(104) <= 0 Then a(112) = 0: GoTo 2925 2920 a(112) = 2.587 * Log(a(104)) 2925 If a(105) <= 0 Then a(113) = 0: GoTo 2935 2930 a(113) = 1.2935 * Log(a(105)) 2935 If a(106) <= 0 Then a(114) = 0: GoTo 2945 2940 a(114) = 1.2935 * Log(a(106)) 2945 If a(107) <= 0 Then a(115) = 0: GoTo 3000 2950 a(115) = 1.2935 * Log(a(107)) 3000 Print #1, "Calcium Oxalate"; Tab(20); Format(a(100), "###0.00"), Tab(43); Format(a(108), "###0.00") 3010 Print #1, "Brushite"; Tab(20); Format(a(101), "###0.00"), Tab(43); Format(a(109), "###0.00") 3020 Print #1, "Hydroxyapatite"; Tab(20); Format(a(102), "###0.00"), Tab(43); Format(a(110), "###0.00") 3030 Rem Print #1, "STRU", Format(A(103), "###0.00"), Tab(43); Format(A(111), "###0.00") 3040 Print #1, "Uric Acid"; Tab(20); Format(a(104), "###0.00"), Tab(43); Format(a(112), "###0.00") 3050 Print #1, "Sodium Urate"; Tab(20); Format(a(105), "###0.00"), Tab(43); Format(a(113), "###0.00") 3060 Print #1, "Ammonium Urate"; Tab(20); Format(a(106), "###0.00"), Tab(43); Format(a(114), "###0.00") 3070 Rem Print #1, "KU", Format(A(100), "###0.00"), Tab(43); Format(A(115), "###0.00") 3080 Print #1,: Print #1,: Print #1,: Print #1, 3100 Close #1 3240 Close 3250 End End Sub Sub Text1_Change() a(1) = Val(Text1.Text) End Sub Sub Text10_Change() a(24) = Val(Text10.Text) End Sub Sub Text11_Change() a(10) = Val(Text11.Text) End Sub Private Sub Text14_Change() p = Text14.Text End Sub Private Sub Text15_Change() q = Text15.Text End Sub Sub Text2_Change() a(2) = Val(Text2.Text) End Sub Sub Text3_Change() a(3) = Val(Text3.Text) End Sub Sub Text4_Change() a(4) = Val(Text4.Text) End Sub Sub Text5_Change() a(91) = Val(Text5.Text) End Sub Sub Text6_Change() a(6) = Val(Text6.Text) End Sub Sub Text7_Change() a(7) = Val(Text7.Text) End Sub Sub Text8_Change() a(9) = Val(Text8.Text) End Sub Sub Text9_Change() a(8) = Val(Text9.Text) End Sub ``` # Note on phosphate and sulfate units (issue #2) Lines `717` and `719` of the V5 source above divide the "Phosphate (mg/dl)" and "Sulfate (mg/dl)" inputs by `3097` and `3206` respectively: ``` 717 a(6) = a(6) / 3097 719 a(7) = a(7) / 3206 ``` These factors are the atomic weights of *inorganic phosphorus* (P, 30.97 g/mol) and *sulfur* (S, 32.06 g/mol) multiplied by 100 to convert mg/dL to mol/L. They do **not** correspond to the molecular weights of the phosphate ion (PO₄³⁻, 94.97 g/mol) or the sulfate ion (SO₄²⁻, 96.07 g/mol). In other words, the V5 algorithm expects "Phosphate (mg/dl)" to be the mass of inorganic phosphorus and "Sulfate (mg/dl)" to be the mass of sulfur — which is consistent with how U.S. clinical labs report those species, but conflicts with the chemical names of the inputs. This mismatch was reported in [issue #2](https://github.com/billdenney/equil2/issues/2) by Lea Lerose, who demonstrated that the R port's unit-aware inputs (e.g. `set_units(32, "mmol_phosphate/L")`) were being inflated by a factor of roughly 94.97 / 30.97 ≈ 3.07 because the `units` package converts the input value via the *true* PO₄ molecular weight while the algorithm then divides by the P-based factor. In the R port, this has been resolved by replacing the V5 factors with the true polyatomic-ion factors: | V5 source line | V5 value | R port value | Reason | |----------------|----------|-------------------|--------------------------------------------------| | `a(6) = a(6) / 3097` | 3097 | **9497** | 94.97 g/mol × 100 dL/L — actual PO₄ ion mass | | `a(7) = a(7) / 3206` | 3206 | **9607** | 96.07 g/mol × 100 dL/L — actual SO₄ ion mass | | output `a(6) * 3097` | 3097 | **9497** | mirror of the input factor | | output `a(7) * 3206` | 3206 | **9607** | mirror of the input factor | Users with legacy clinical data recorded as "phosphate (mg/dL)" meaning *mass of P* should convert by multiplying by ≈ 3.066 (= 94.97 / 30.97) to get mass of PO₄, or — more simply — supply the value in mmol/L, where 1 mmol of P = 1 mmol of PO₄. Same logic applies to sulfate with the ratio 96.07 / 32.06 ≈ 3.00. **Cross-validation.** The R port has been cross-validated against the V5 BASIC source itself by transcribing the `calc()` Sub above (with the same `/9497` and `/9607` fix applied) into LibreOffice 26.2.4 Basic and running it with the LabCorp inputs in `mmol_phosphate/L` / `mmol_sulfate/L` form. R and V5 BASIC agree on ionic strength to about 2 parts per billion and on all six supersaturation values (Calcium Oxalate, Brushite, Hydroxyapatite, Uric Acid, Sodium Urate, Ammonium Urate) plus their ΔGibbs energies to better than 1 part per million. Residual differences are floating-point ordering noise. # V1 vs V5 differences (informational) The R port is based on the V5 (form-based BASIC) source shown above. The original Mayo Clinic class module `clsEquil2` ("V1") is a richer implementation that is provided in issue #2 as an attached `.txt` file. Below is a brief summary of the differences, included so users can decide whether they need V1's additional capabilities. **Inputs.** V1 takes all species in mmol per collection volume (plus a `Vol` parameter) and then converts to mol/L by dividing by `Vol × 1000`. This means V1 has no per-species atomic-weight factors at all, and the phosphate/sulfate confusion above does not arise — the V1 approach is unambiguous by construction. **Additional chemistry in V1.** V1 models several pathways that V5 (and therefore the R port) does not: - **CO₂ / bicarbonate / carbonate** chemistry (H₂CO₃, HCO₃⁻, CO₃²⁻, plus NaCO₃, Na₂CO₃, CaCO₃, MgCO₃ complexes). - **Pyrophosphate (PP)** speciation (HPP, H₂PP, H₃PP, H₄PP, CaPP, CaHPP, CaOHPP, MgPP, MgOHPP, NaPP, Na₂PP, NaHPP, KPP). - **TRIS buffer** (HTRIS). - **Struvite (MgNH₄PO₄)** included in the Mg and NH₄ mass balances. - **Na-urate / NH₄-urate / K-urate** complexes included in the Na/NH₄/K mass balances (V5 computes their supersaturation ratios but does not feed them back into the cation balance). - **Chloride auto-fill** via electroneutrality when no chloride is supplied (V5 requires a chloride input). **Stability constants.** V1 and V5 share the same constants for the core CaOx / brushite / hydroxyapatite / uric-acid / Na-urate / NH₄-urate pathways. The constants that differ are concentrated in the CO₂, PP, and HU pathways — pathways V5 does not use anyway. For example, V1 uses `P7 = 0.00229` for CO₂ dissociation while V5 uses `S(61) = 0.00123`; V1 uses `7.943 × 10¹⁰` for HU formation while V5 uses `S(9) = 278000`. **Outputs.** V1 reports activity products (`APCaOx`, `APBr`, `APStru`, `APUA`, `APNaU`) and relative saturation ratios (RSR) for calcium oxalate, brushite, struvite, and uric acid only. V5 (and the R port) additionally reports supersaturation for hydroxyapatite, sodium urate, and ammonium urate, plus negative ΔGibbs energies for all six. **Convergence.** V1 caps iteration at 500 and treats non-convergence as an error. V5 (and the R port) caps at 50 and warns. **Bottom line.** The R port is a faithful translation of V5. The only true *bug* it inherited from V5 is the phosphate/sulfate unit confusion fixed in this release. The other items above are V5-vs-V1 *scope* gaps, not translation errors. A more complete port that adopts V1's structure is feasible but is out of scope for the current release.