VERSION 5.00
Begin VB.Form Form1 
   BackColor       =   &H8000000A&
   Caption         =   "Energy-gap Law Curves"
   ClientHeight    =   7680
   ClientLeft      =   45
   ClientTop       =   330
   ClientWidth     =   7320
   LinkTopic       =   "Form1"
   ScaleHeight     =   7680
   ScaleWidth      =   7320
   StartUpPosition =   3  'Windows Default
   Begin VB.TextBox Text11 
      Height          =   285
      Left            =   6360
      TabIndex        =   41
      Text            =   "298"
      Top             =   1560
      Width           =   855
   End
   Begin VB.HScrollBar HScroll8 
      Height          =   255
      LargeChange     =   10
      Left            =   5400
      Max             =   320
      Min             =   260
      TabIndex        =   39
      Top             =   1320
      Value           =   298
      Width           =   1815
   End
   Begin VB.OptionButton Option4 
      Caption         =   "Marcus-Bronsted"
      Height          =   375
      Left            =   3840
      TabIndex        =   38
      Top             =   1440
      Width           =   1575
   End
   Begin VB.OptionButton Option3 
      Caption         =   "Marcus"
      Height          =   375
      Left            =   2880
      TabIndex        =   37
      Top             =   1440
      Width           =   855
   End
   Begin VB.OptionButton Option2 
      Caption         =   "MD-Bronsted"
      Height          =   375
      Left            =   1440
      TabIndex        =   36
      Top             =   1440
      Width           =   1335
   End
   Begin VB.OptionButton Option1 
      Caption         =   "Moser-Dutton"
      Height          =   375
      Left            =   0
      TabIndex        =   35
      Top             =   1440
      Width           =   1455
   End
   Begin VB.TextBox Text10 
      Height          =   288
      Left            =   6720
      TabIndex        =   33
      Top             =   240
      Width           =   492
   End
   Begin VB.TextBox Text9 
      Height          =   288
      Left            =   2880
      TabIndex        =   32
      Text            =   "Text9"
      Top             =   720
      Width           =   492
   End
   Begin VB.TextBox Text8 
      Height          =   288
      Left            =   2880
      TabIndex        =   31
      Text            =   "Text8"
      Top             =   360
      Width           =   492
   End
   Begin VB.CommandButton Command5 
      Caption         =   "Save"
      Height          =   252
      Left            =   6600
      TabIndex        =   30
      Top             =   840
      Width           =   612
   End
   Begin VB.ComboBox Combo1 
      Height          =   288
      Left            =   2280
      TabIndex        =   28
      Text            =   "Combo1"
      Top             =   1080
      Width           =   2172
   End
   Begin VB.CommandButton Command4 
      Caption         =   "Zoom"
      Height          =   252
      Left            =   6600
      TabIndex        =   27
      Top             =   1080
      Width           =   612
   End
   Begin VB.CommandButton Command3 
      Caption         =   "Print"
      Height          =   252
      Left            =   5760
      TabIndex        =   26
      Top             =   1080
      Width           =   612
   End
   Begin VB.HScrollBar HScroll6 
      Height          =   252
      LargeChange     =   500
      Left            =   1080
      Max             =   3000
      Min             =   10
      SmallChange     =   10
      TabIndex        =   25
      Top             =   720
      Value           =   500
      Width           =   1332
   End
   Begin VB.TextBox Text7 
      Height          =   288
      Left            =   5760
      TabIndex        =   22
      Text            =   "6.3"
      Top             =   720
      Width           =   492
   End
   Begin VB.HScrollBar HScroll7 
      Height          =   252
      LargeChange     =   100
      Left            =   4200
      Max             =   3000
      SmallChange     =   5
      TabIndex        =   21
      Top             =   720
      Value           =   630
      Width           =   1452
   End
   Begin VB.TextBox Text6 
      Height          =   288
      Left            =   2400
      TabIndex        =   19
      Text            =   "0.5"
      Top             =   720
      Width           =   492
   End
   Begin VB.CommandButton Command2 
      Caption         =   "Param vals"
      Height          =   252
      Left            =   4560
      TabIndex        =   17
      Top             =   1080
      Width           =   1092
   End
   Begin VB.TextBox Text5 
      Height          =   288
      Left            =   6720
      TabIndex        =   16
      Text            =   "-275"
      Top             =   0
      Width           =   492
   End
   Begin VB.HScrollBar HScroll5 
      Height          =   252
      LargeChange     =   10
      Left            =   5520
      Max             =   100
      Min             =   -650
      TabIndex        =   15
      Top             =   0
      Value           =   -275
      Width           =   1092
   End
   Begin VB.TextBox Text4 
      Height          =   288
      Left            =   4200
      TabIndex        =   13
      Text            =   "-90"
      Top             =   0
      Width           =   372
   End
   Begin VB.TextBox Text3 
      Height          =   288
      Left            =   1920
      TabIndex        =   12
      Text            =   "300"
      Top             =   0
      Width           =   372
   End
   Begin VB.HScrollBar HScroll4 
      Height          =   252
      LargeChange     =   10
      Left            =   3000
      Max             =   0
      Min             =   -250
      TabIndex        =   9
      Top             =   0
      Value           =   -90
      Width           =   1212
   End
   Begin VB.HScrollBar HScroll3 
      Height          =   252
      LargeChange     =   10
      Left            =   720
      Max             =   600
      Min             =   1
      TabIndex        =   8
      Top             =   0
      Value           =   300
      Width           =   1212
   End
   Begin VB.CommandButton Command1 
      Caption         =   "Clear"
      Height          =   252
      Left            =   6600
      TabIndex        =   7
      Top             =   600
      Width           =   612
   End
   Begin VB.TextBox Text2 
      Height          =   288
      Left            =   5760
      TabIndex        =   4
      Text            =   "2.6"
      Top             =   360
      Width           =   492
   End
   Begin VB.TextBox Text1 
      Height          =   288
      Left            =   2400
      TabIndex        =   3
      Text            =   "2.0"
      Top             =   360
      Width           =   492
   End
   Begin VB.HScrollBar HScroll2 
      Height          =   252
      LargeChange     =   100
      Left            =   4200
      Max             =   3000
      SmallChange     =   5
      TabIndex        =   2
      Top             =   360
      Value           =   260
      Width           =   1452
   End
   Begin VB.HScrollBar HScroll1 
      Height          =   252
      LargeChange     =   500
      Left            =   1080
      Max             =   3500
      Min             =   10
      SmallChange     =   10
      TabIndex        =   1
      Top             =   360
      Value           =   2630
      Width           =   1332
   End
   Begin VB.PictureBox Picture1 
      AutoRedraw      =   -1  'True
      BackColor       =   &H8000000A&
      BeginProperty Font 
         Name            =   "Times New Roman"
         Size            =   10.5
         Charset         =   0
         Weight          =   400
         Underline       =   0   'False
         Italic          =   0   'False
         Strikethrough   =   0   'False
      EndProperty
      Height          =   5655
      Left            =   0
      ScaleHeight     =   5595
      ScaleWidth      =   7155
      TabIndex        =   0
      Top             =   1920
      Width           =   7215
   End
   Begin VB.Label Label12 
      Caption         =   "Temp K"
      Height          =   255
      Left            =   5400
      TabIndex        =   40
      Top             =   1560
      Width           =   735
   End
   Begin VB.Label Label11 
      Caption         =   "  logKs"
      Height          =   252
      Left            =   6240
      TabIndex        =   34
      Top             =   240
      Width           =   492
   End
   Begin VB.Label Label10 
      Caption         =   "Change hidden parameters"
      Height          =   252
      Left            =   120
      TabIndex        =   29
      Top             =   1080
      Width           =   2052
   End
   Begin VB.Label Label9 
      BackColor       =   &H0000FFFF&
      Caption         =   "2nd e-"
      Height          =   252
      Left            =   0
      TabIndex        =   24
      Top             =   720
      Width           =   492
   End
   Begin VB.Label Label8 
      BackColor       =   &H000000FF&
      Caption         =   "1st e-"
      Height          =   252
      Left            =   0
      TabIndex        =   23
      Top             =   360
      Width           =   492
   End
   Begin VB.Label Label7 
      Caption         =   "distance"
      Height          =   252
      Left            =   3480
      TabIndex        =   20
      Top             =   720
      Width           =   732
   End
   Begin VB.Label Label6 
      Caption         =   "  lamda"
      Height          =   252
      Left            =   480
      TabIndex        =   18
      Top             =   720
      Width           =   612
   End
   Begin VB.Label Label5 
      Caption         =   "  delta Em"
      Height          =   252
      Left            =   4680
      TabIndex        =   14
      Top             =   0
      Width           =   732
   End
   Begin VB.Label Label4 
      Caption         =   "Em bL"
      Height          =   252
      Left            =   2400
      TabIndex        =   11
      Top             =   0
      Width           =   492
   End
   Begin VB.Label Label3 
      Caption         =   "Em ISP"
      Height          =   252
      Left            =   0
      TabIndex        =   10
      Top             =   0
      Width           =   612
   End
   Begin VB.Label Label2 
      Caption         =   "distance"
      Height          =   252
      Left            =   3480
      TabIndex        =   6
      Top             =   360
      Width           =   612
   End
   Begin VB.Label Label1 
      Caption         =   "  lamda"
      Height          =   252
      Left            =   480
      TabIndex        =   5
      Top             =   360
      Width           =   612
   End
End
Attribute VB_Name = "Form1"
Attribute VB_GlobalNameSpace = False
Attribute VB_Creatable = False
Attribute VB_PredeclaredId = True
Attribute VB_Exposed = False
Dim Newline$, zoomtog%, Occup, Fix_Ks%, Fix_deltaE%, OldEmISP, logKs
Dim testEact1, testEact2, Lamda1cutoff, Lamda2cutoff, Gact
Dim lamda1, dist1, lamda2, dist2, EmISP, EmbL, EmQQH2, EmQQH2bound
Dim deltaG, deltaGprot, deltaE1elec
Dim deltaE, deltaE1, deltaE2, deltaE2mV, deltaG_bound, slope_logk1, dlogk1
Dim EmQHdotQH2, EmQQHdot, deltaG_overall, Rate1st, Rate2nd, lastval1
Dim pKISPox, pKQH2, beta, Rgc, Farad, TK, MD, MDB, Marcus, MB, Fac, RTV, RTmV

Private Sub Combo1_Click()
' Change default values for hidden parameters
If Combo1.ListIndex = 0 Then
    EmQQH2 = Val(InputBox("Change Em(Q/QH2)?", "Current value is " & Str$(EmQQH2), EmQQH2))
    Combo1.List(0) = "Em(Q/QH2) = " & Str$(EmQQH2)
ElseIf Combo1.ListIndex = 1 Then
    EmQQH2bound = Val(InputBox("Change Em(Q/QH2)(bound)?", "Current value is " & Str$(EmQQH2bound), EmQQH2bound))
    Combo1.List(1) = "Em(Q/QH2)(bound) = " & Str$(EmQQH2bound)
ElseIf Combo1.ListIndex = 2 Then
    testEact1 = Val(InputBox("Change Eact 1st e?", "Current value is " & Str$(testEact1), testEact1))
    Combo1.List(2) = "Eact 1st e = " & Str$(testEact1)
    CalcLamda1Cutoff
ElseIf Combo1.ListIndex = 3 Then
    testEact2 = Val(InputBox("Change Eact 2nd e?", "Current value is " & Str$(testEact2), testEact2))
    Combo1.List(3) = "Eact 2nd e = " & Str$(testEact2)
    CalcLamda2Cutoff
ElseIf Combo1.ListIndex = 4 Then
    pKISPox = Val(InputBox("Change pK_ISPox?", "Current value is " & Str$(pKISPox), pKISPox))
    Combo1.List(4) = "pKISPox = " & Str$(pKISPox)
ElseIf Combo1.ListIndex = 5 Then
    pKQH2 = Val(InputBox("Change pK_QH2?", "Current value is " & Str$(pKQH2), pKQH2))
    Combo1.List(5) = "pKQH2 = " & Str$(pKQH2)
End If
CalcGvals
End Sub

Private Sub Command1_Click()
    Picture1.Cls
End Sub
Static Function Log10(X)
    Log10 = Log(X) / Log(10#)
End Function

Private Sub Command2_Click()
'Print Parameters
Msg0 = "Em(ISP) = " & Str$(EmISP) & " mV; Em(bL) = " & Str$(EmbL) & " mV; delta E, 1st e- = " & Str$(deltaE) & " mV; delta E, 2nd e- = " & Str$(deltaE2mV)
Msg1 = "EmQ/QH2 (pool) = " & Str$(EmQQH2) & ", EmQ/QH2 (bound) = " & Str$(EmQQH2bound) & ", EmQ/QH. = " & Str$(EmQQHdot) & ", EmQH./QH2 = " & Str$(EmQHdotQH2)
Msg2 = "lamda, 1st = " & Str$(lamda1) & ", distance 1st = " & Str$(dist1) & ", lamda, 2nd = " & Str$(lamda2) & "; distance 2nd = " & Str$(dist2) & "; Occupancy = " & Format(Occup, "#.##e-##")
Msg3 = "Rate, 1st e- = " & Format(Rate1st, "#.##e-##") & " /s; Rate 2nd e- = " & Format(Rate2nd, "#.##e-##") & "; Slope at E1 = " & Format(slope_logk1, "#.####") & " /mV" & ", pKISPox = " & Str$(pKISPox) & ", pKQH2 = " & Str$(pKQH2)

MsgBox (Msg0 + Newline$ + Msg1 + Newline$ + Msg2 + Newline$ + Msg3)
End Sub

Private Sub Command3_Click()
'Print Graph
Msg0 = "Em(ISP) = " & Str$(EmISP) & " mV; Em(bL) = " & Str$(EmbL) & " mV"
Msg05 = "delta E, 1st e- = " & Str$(deltaE) & " mV; delta E, 2nd e- = " & Str$(deltaE2mV)
Msg1 = "EmQ/QH2 (pool) = " & Str$(EmQQH2) & ", EmQ/QH2 (bound) = " & Str$(EmQQH2bound)
Msg15 = "EmQ/QH. = " & Str$(EmQQHdot) & ", EmQH./QH2 = " & Str$(EmQHdotQH2) & "; Occupancy = " & Format(Occup, "#.##e-##")
Msg2 = "lamda, 1st = " & Str$(lamda1) & ", distance 1st = " & Str$(dist1)
Msg25 = "lamda, 2nd = " & Str$(lamda2) & "; distance 2nd = " & Str$(dist2)
Msg3 = "Rate 1st e = " & Format(Rate1st, "#.##e-##") & "/s; Rate 2nd e = " & Format(Rate2nd, "#.##e-##") & "/s)"
Msg35 = "Slope (E1) = " & Format(slope_logk1, "#.####") & "/mV"
If MDB = True Or MB = True Then
    Msg35 = Msg35 & ", pK, ISPox = " & Str$(pKISPox) & "; pKQH2 = " & Str$(pKQH2) & ""
End If


Printer.Print
Printer.Print
Printer.Print "Parameters for this simulation"
Printer.Print
Printer.Print Msg0
Printer.Print
Printer.Print Msg05
Printer.Print
Printer.Print Msg1
Printer.Print
Printer.Print Msg15
Printer.Print
Printer.Print Msg2
Printer.Print
Printer.Print Msg25
Printer.Print
Printer.Print Msg3
Printer.Print
Printer.Print Msg35

Printer.NewPage
Printer.EndDoc
If zoomtog% = 0 Then
    Printer.Scale (3, 25)-(-4, -17.5)
    Printer.Line (2, 15)-(-3, -5), , B
ElseIf zoomtog% = 1 Then
    Printer.Scale (2, 24)-(-3, -12.5)
    Printer.Line (1, 14)-(-2, 0), , B
End If

CalcGvals

'x-axis tick marks
If zoomtog% = 0 Then
    For tmp = -3 To 2
        Printer.Line (tmp, -5)-(tmp, -4.5)
        Printer.PSet (tmp, -5.5), BackColor
        Printer.Print tmp
    Next
    For tmp = -3 To 2 Step 0.2
        Printer.Line (tmp, -5)-(tmp, -4.8)
    Next
ElseIf zoomtog% = 1 Then
    For tmp = -2 To 1
        Printer.Line (tmp, 0)-(tmp, 0.5)
        Printer.PSet (tmp, -0.5), BackColor
        Printer.Print tmp
    Next
    For tmp = -2 To 1 Step 0.2
        Printer.Line (tmp, 0)-(tmp, 0.2)
    Next
End If

'Print x-axis label
If zoomtog% = 0 Then
    Printer.PSet (-0.1, -6.5), BackColor
    Printer.Print "delta G (V)"
ElseIf zoomtog% = 1 Then
    Printer.PSet (-0.1, -1), BackColor
    Printer.Print "delta G (V)"
End If
'y-axis tick marks
If zoomtog% = 0 Then
    For tmp = -5 To 15
        If tmp Mod 5 Then   'minor tics
            Printer.Line (-3, tmp)-(-2.95, tmp)
            Printer.Line (2, tmp)-(1.95, tmp)
        Else
            'major tics
            'right axis
            Printer.Line (-3, tmp)-(-2.875, tmp)
            Printer.PSet (-3.1, tmp), BackColor
            out$ = Format(10 ^ tmp, "#.0 E-##")
            Printer.Print out$
            'left axis
            Printer.Line (2, tmp)-(1.875, tmp)
            Printer.PSet (2.3, tmp), BackColor
            Printer.Print tmp
        End If
    Next
ElseIf zoomtog% = 1 Then
    For tmp = 0 To 15
        If tmp Mod 5 Then   'minor tics
            Printer.Line (-2, tmp)-(-1.95, tmp)
            Printer.Line (1, tmp)-(0.95, tmp)
        Else
            'major tics
            'right axis
            Printer.Line (-2, tmp)-(-1.875, tmp)
            Printer.PSet (-2.1, tmp), BackColor
            out$ = Format(10 ^ tmp, "#.0 E-##")
            Printer.Print out$
            'left axis
            Printer.Line (1, tmp)-(0.875, tmp)
            Printer.PSet (1.3, tmp), BackColor
            Printer.Print tmp
        End If
    Next
End If
'Print y-axis labels
If zoomtog% = 0 Then
    Printer.PSet (-3.5, 7), BackColor
    Printer.Print "k (/s)"
    Printer.PSet (2.85, 7), BackColor
    Printer.Print "log k"
ElseIf zoomtog% = 1 Then
    Printer.PSet (-2.5, 7), BackColor
    Printer.Print "k (/s)"
    Printer.PSet (1.85, 7), BackColor
    Printer.Print "log k"
End If

tg% = 0
ddG = 0
lastdG = 0
lastval1 = 0
lastval2 = 0
overslope = 0

'Calculate Occupancy of semquinone
K1 = 10 ^ (deltaE / RTmV)
Occup = K1 / (1 + K1)

'Marcus curves
'Ist electron
For deltaG = -3.5 * lamda1 To 1.5 * lamda1 Step 0.01  '10 mV steps
    If MD = True Then
        logk1 = 15 - 0.6 * dist1 - (3.1 * 298 / TK) * ((deltaG + lamda1) ^ 2 / lamda1)
    ElseIf MDB = True Then
        logk1 = 15 - 0.6 * dist1 - (3.1 * 298 / TK) * ((deltaG + lamda1) ^ 2 / lamda1) - (pKQH2 - pKISPox)
    ElseIf Marcus = True Then
        logk1 = 13 - (beta * (dist1 - 3.6)) / Log(10) - (Fac / TK) * (deltaG + lamda1) ^ 2 / lamda1
    ElseIf MB = True Then
        logk1 = 13 - (beta * (dist1 - 3.6)) / Log(10) - ((Fac / TK) * (deltaG + lamda1) ^ 2 / lamda1) - (pKQH2 - pKISPox)
    End If
'    logk1 = 15 - 0.6 * dist1 - 3.1 * (deltaG + lamda1) ^ 2 / lamda1
        If tg% = 0 Then
            If logk1 >= -5 Then
                If deltaG >= -3 Then
                    Printer.PSet (deltaG, logk1)
                    tg% = 1
                    lastval1 = logk1
                    lastdG = deltaG
                End If
            End If
        ElseIf logk1 >= -5 Then
            If deltaG >= -3 Then
                EtpRate = 10 ^ logk1
                    If EtpRate > 1.5 * 10 ^ 3 Then
                        Printer.DrawStyle = 1
                    Else
                        Printer.DrawStyle = 0
                    End If
                Printer.Line -(deltaG, logk1), LineCol
                'Find change in logk for this 10 mV step
                dlogk1 = logk1 - lastval1
                ddG = deltaG - lastdG
                'Look for intercept with deltaE values for 1st electron transfer
                If MD = True Or Marcus = True Then
                    If (-deltaE1 > lastdG And -deltaE1 < deltaG) Or -deltaE1 = deltaG Then
                        slope_logk1 = dlogk1 / 10
                        Printer.PSet (deltaG, logk1)
                        Printer.Print Format(10 ^ logk1, "#.##e-##"); "/s, "; Format(slope_logk1, "##.####") & "/mV"
                        Printer.PSet (deltaG, logk1)
                        Rate1st = 10 ^ logk1
                    End If
                ElseIf MDB = True Or MB = True Then
                   If (-deltaE1elec > lastdG And -deltaE1elec < deltaG) Or -deltaE1elec = deltaG Then
                        slope_logk1 = dlogk1 / 10
                        Printer.PSet (deltaG, logk1)
                        Printer.Print Format(10 ^ logk1, "#.##e-##"); "/s, "; Format(slope_logk1, "##.####") & "/mV"
                        Printer.PSet (deltaG, logk1)
                        Rate1st = 10 ^ logk1
                    End If
                End If
                lastval1 = logk1
                lastdG = deltaG
            End If
        End If
Next
tg% = 0
overslope = 0
'Second electron
For deltaG = -3.5 * lamda2 To 1.5 * lamda2 Step 0.01  '10 mV steps
    If MD = True Then
        logk2 = 15 - 0.6 * dist2 - (3.1 * 298 / TK) * (deltaG + lamda2) ^ 2 / lamda2
    ElseIf MDB = True Then
        logk2 = 15 - 0.6 * dist2 - (3.1 * 298 / TK) * (deltaG + lamda2) ^ 2 / lamda2
    ElseIf Marcus = True Then
        logk2 = 13 - (beta * (dist2 - 3.6)) / Log(10) - (Fac / TK) * (deltaG + lamda2) ^ 2 / lamda2
    ElseIf MB = True Then
        logk2 = 13 - (beta * (dist2 - 3.6)) / Log(10) - (Fac / TK) * (deltaG + lamda2) ^ 2 / lamda2
    End If
'    logk2 = 15 - 0.6 * dist2 - 3.1 * (deltaG + lamda2) ^ 2 / lamda2
        If tg% = 0 Then
            If logk2 >= -5 Then
                If deltaG >= -3 Then
                    Printer.PSet (deltaG, logk2)
                    tg% = 1
                    lastval2 = logk2
                    lastdG = deltaG
                End If
            End If
        ElseIf logk2 >= -5 Then
            If deltaG >= -3 Then
                EtpRate = Occup * 10 ^ logk2
                    If EtpRate > 1.5 * 10 ^ 3 Then
                        Printer.DrawStyle = 1
                    Else
                        Printer.DrawStyle = 0
                    End If
                Printer.Line -(deltaG, logk2), LineCol
                If (-deltaE2 > lastdG And -deltaE2 < deltaG) Or -deltaE2 = deltaG Then
                    Printer.PSet (deltaG, logk2)
                    Printer.Print Format(Occup * 10 ^ logk2, "#.##E-##")
                    Printer.PSet (deltaG, logk2)
                    Rate2nd = Occup * 10 ^ logk2
                End If
                dlogk2 = logk2 - lastval2
                lastval2 = logk2
                lastdG = deltaG
            End If
        End If
Next

If MD = True Or Marcus = True Then
    Printer.Line (-deltaE1, -3)-(-deltaE1, 12.5)
    Printer.PSet (-deltaE1, 13), BackColor
    Printer.Print "E1"
ElseIf MDB = True Or MB = Trur Then
    Printer.Line (-deltaE1elec, -3)-(-deltaE1elec, 12.5)
    Printer.PSet (-deltaE1elec, 13), BackColor
    Printer.Print "E1'"
End If
Printer.Line (-deltaE2, -3)-(-deltaE2, 13)
Printer.PSet (-deltaE2, 13.5), BackColor
Printer.Print "E2"
Printer.Line (deltaG_bound, -3)-(deltaG_bound, 13.5)
Printer.PSet (deltaG_bound, 14), BackColor
Printer.Print "Gb"
Printer.Line (deltaG_overall, -3)-(deltaG_overall, 14)
Printer.PSet (deltaG_overall, 14.5), BackColor
Printer.Print "Go"
Printer.DrawStyle = 1
Printer.Line (-3, Log10(1500))-(2, Log10(1500))
Printer.DrawStyle = 0
Printer.EndDoc
End Sub

Private Sub Command4_Click()
If zoomtog% = 0 Then
    zoomtog% = 1
Else
    zoomtog% = 0
End If
End Sub

Private Sub Command5_Click()
'Save
Msg0 = "Em(ISP) = " & Str$(EmISP) & " mV; Em(bL) = " & Str$(EmbL) & " mV; deltaGe, 1st e- = " & Str$(deltaE) & " mV; deltaGe, 2nd e- = " & Str$(deltaE2mV)
Msg1 = "EmQ/QH2 (pool) = " & Str$(EmQQH2) & ", EmQ/QH2 (bound) = " & Str$(EmQQH2bound) & ", EmQ/QH. = " & Str$(EmQQHdot) & ", EmQH./QH2 = " & Str$(EmQHdotQH2)
Msg2 = "lamda, 1st = " & Str$(lamda1) & ", distance 1st = " & Str$(dist1) & ", lamda, 2nd = " & Str$(lamda2) & "; distance 2nd = " & Str$(dist2)
Msg3 = "Rate, 1st e- = " & Format(Rate1st, "#.##e-##") & " /s; Rate 2nd e- = " & Format(Rate2nd, "#.##e-##") & "; Slope at E1 (or E1') = " & Format(slope_logk1, "#.####") & " /mV"
If MDB = True Or MB = True Then
    Msg4 = "pK, ISPox = " & Str$(pKISPox) & "; pKQH2 = " & Str$(pKQH2) & "; deltaE 1st e-transfer = " & Str$(deltaE1elec * 1000) & " mV"
Else: Msg4 = ""
End If

Filenam$ = InputBox("Name for file to save this set?", "File data", "Marcus_data.txt")
    Open Filenam$ For Output As #1
    Write #1, Msg0, Msg1, Msg2, Msg3, Msg4
    
    'Marcus curves

'1st electron

Write #1,
Write #1, "1st electron"
Write #1,

For deltaG = -3.5 * lamda1 To 1.5 * lamda1 Step 0.01  '10 mV steps
    
    If MD = True Then
        logk1 = 15 - 0.6 * dist1 - (3.1 * 298 / TK) * ((deltaG + lamda1) ^ 2 / lamda1)
    ElseIf MDB = True Then
        logk1 = 15 - 0.6 * dist1 - (3.1 * 298 / TK) * ((deltaG + lamda1) ^ 2 / lamda1) - (pKQH2 - pKISPox)
    ElseIf Marcus = True Then
        logk1 = 13 - (beta * (dist1 - 3.6)) / Log(10) - (Fac / TK) * (deltaG + lamda1) ^ 2 / lamda1
    ElseIf MB = True Then
        logk1 = 13 - (beta * (dist1 - 3.6)) / Log(10) - ((Fac / TK) * (deltaG + lamda1) ^ 2 / lamda1) - (pKQH2 - pKISPox)
    End If

        If tg% = 0 Then
            If logk1 >= -5 Then
                If deltaG >= -3 Then
                    tg% = 1
                    lastval1 = logk1
                    lastdG = deltaG
                End If
            End If
        ElseIf logk1 >= -5 Then
            If deltaG >= -3 Then
                Write #1, Int(deltaG * 1000) / 1000, Int(logk1 * 1000) / 1000
            End If
        End If
Next
tg% = 0
overslope = 0
'Second electron
Write #1,
Write #1, "2nd electron"
Write #1,

For deltaG = -3.5 * lamda2 To 1.5 * lamda2 Step 0.01  '10 mV steps
    'logk2 = 15 - 0.6 * dist2 - 3.1 * (deltaG + lamda2) ^ 2 / lamda2
    
    If MD = True Then
        logk2 = 15 - 0.6 * dist2 - (3.1 * 298 / TK) * (deltaG + lamda2) ^ 2 / lamda2
    ElseIf MDB = True Then
        logk2 = 15 - 0.6 * dist2 - (3.1 * 298 / TK) * (deltaG + lamda2) ^ 2 / lamda2
    ElseIf Marcus = True Then
        logk2 = 13 - (beta * (dist2 - 3.6)) / Log(10) - (Fac / TK) * (deltaG + lamda2) ^ 2 / lamda2
    ElseIf MB = True Then
        logk2 = 13 - (beta * (dist2 - 3.6)) / Log(10) - (Fac / TK) * (deltaG + lamda2) ^ 2 / lamda2
    End If
    
        If tg% = 0 Then
            If logk2 >= -5 Then
                If deltaG >= -3 Then
                     tg% = 1
                    lastval2 = logk2
                    lastdG = deltaG
                End If
            End If
        ElseIf logk2 >= -5 Then
            If deltaG >= -3 Then
                Write #1, Int(deltaG * 1000) / 1000, Int(logk2 * 1000) / 1000
            End If
        End If
Next
Close #1
End Sub

Private Sub Form_Load()
'Set up starting parameters

Newline$ = Chr$(13) & Chr$(10)

Combo1.AddItem "Em(Q/QH2) = 90", 0
Combo1.AddItem "Em(Q/QH2)(bound) = 130", 1
Combo1.AddItem "Eact 1st e = 0.65", 2
Combo1.AddItem "Eact 2nd e = 0.26", 3
Combo1.AddItem "pKISPox = 7.6", 4
Combo1.AddItem "pKQH2 = 11.5", 5
Combo1.Text = "Default values"
EmQQH2 = 90
EmQQH2bound = 130

HScroll1.Value = 2010
lamda1 = 2.01
Text1.Text = 2.01

HScroll2.Value = 600
dist1 = 6#
Text2.Text = 6#

HScroll6.Value = 1200
lamda2 = 1.2
Text6.Text = 1.2

HScroll7.Value = 630
dist2 = 6.3
Text7.Text = 6.3

Rgc = 8.31447                   'J.K-1.mol-1, Gas Constant, R
Farad = 96485.3383              'Faraday constant, C/mol or J/V
TK = 298                        'T in degrees K
'Fac = 4.22800 * 298 at 298 K
Fac = Farad / (4 * Log(10) * Rgc)    'F/4 * 2.303 *R
RTV = 2.303 * Rgc * TK / Farad  '2.303RT/F in volts
RTmV = RTV * 1000               '2.303RT/F in millivolts
HScroll8.Value = TK
Text11.Text = 298


HScroll3.Value = 310
Text3.Text = 310
EmISP = 310

HScroll4.Value = -90
Text4.Text = -90
EmbL = -90

HScroll5.Value = -275
Text5.Text = -275
deltaE = -275

pKQH2 = 11.5
pKISPox = 7.6
testEact1 = 0.65
CalcLamda1Cutoff
testEact2 = 0.26
CalcLamda2Cutoff

MD = False
MDB = True
Marcus = False
MB = False
Option2.Value = True
beta = 1.4          'slope of Dutton's ruler

Fix_Ks% = 0
Fix_deltaE% = 1
Text10.BackColor = QBColor(15)
Text5.BackColor = QBColor(10)

End Sub

Private Sub HScroll1_Change()
lamda1 = HScroll1.Value / 1000
Text1.Text = lamda1
If lamda1 > Lamda1cutoff Then
    Text1.BackColor = QBColor(11)
Else
    Text1.BackColor = QBColor(15)
End If
End Sub

Private Sub HScroll2_Change()
dist1 = HScroll2.Value / 100
Text2.Text = dist1
End Sub

Private Sub HScroll3_Change()
EmISP = HScroll3.Value
Text3.Text = EmISP
If Fix_Ks% = 1 Then
    Change = OldEmISP - EmISP
    deltaE = deltaE - Change
    HScroll5.Value = deltaE
    Text5.Text = deltaE
ElseIf Fix_deltaE% = 1 Then
    'We leave deltaE the same
End If

CalcGvals
OldEmISP = EmISP
End Sub

Private Sub HScroll4_Change()
EmbL = HScroll4.Value
Text4.Text = EmbL
CalcGvals
End Sub

Private Sub HScroll5_Change()
deltaE = HScroll5.Value
Text5.Text = deltaE
CalcGvals
CalcLamda1Cutoff
If lamda1 > Lamda1cutoff Then
    Text1.BackColor = QBColor(11)
Else
    Text1.BackColor = QBColor(15)
End If
CalcLamda2Cutoff
If lamda2 > Lamda2cutoff Then
    Text6.BackColor = QBColor(11)
Else
    Text6.BackColor = QBColor(15)
End If
End Sub

Private Sub HScroll6_Change()
lamda2 = HScroll6.Value / 1000
Text6.Text = lamda2
If lamda2 > Lamda2cutoff Then
    Text6.BackColor = QBColor(11)
Else
    Text6.BackColor = QBColor(15)
End If
End Sub

Private Sub HScroll7_Change()
dist2 = HScroll7.Value / 100
Text7.Text = dist2

End Sub

Private Sub HScroll8_Change()
    TK = HScroll8.Value
    Text11.Text = TK
End Sub

Private Sub Option1_Click()
    MD = True
    MDB = False
    Marcus = False
    MB = False
End Sub

Private Sub Option2_Click()
    MD = False
    MDB = True
    Marcus = False
    MB = False
End Sub

Private Sub Option3_Click()
    MD = False
    MDB = False
    Marcus = True
    MB = False
End Sub

Private Sub Option4_Click()
    MD = False
    MDB = False
    Marcus = False
    MB = True
End Sub

Private Sub Picture1_Click()

'Draw graph on screen

If zoomtog% = 0 Then
    Picture1.Scale (3, 16.5)-(-4, -7.5)
    Picture1.Line (2, 15)-(-3, -5), , B
ElseIf zoomtog% = 1 Then
    Picture1.Scale (2, 14.5)-(-3, -1.5)
    Picture1.Line (1, 14)-(-2, 0), , B
End If

CalcGvals

'x-axis tick marks
If zoomtog% = 0 Then
    For tmp = -3 To 2
        Picture1.Line (tmp, -5)-(tmp, -4.5)
        Picture1.PSet (tmp, -5.5), BackColor
        Picture1.Print tmp
    Next
    For tmp = -3 To 2 Step 0.2
        Picture1.Line (tmp, -5)-(tmp, -4.8)
    Next
ElseIf zoomtog% = 1 Then
    For tmp = -2 To 1
        Picture1.Line (tmp, 0)-(tmp, 0.5)
        Picture1.PSet (tmp, -0.5), BackColor
        Picture1.Print tmp
    Next
    For tmp = -2 To 1 Step 0.2
        Picture1.Line (tmp, 0)-(tmp, 0.2)
    Next
End If

'Print x-axis label
If zoomtog% = 0 Then
    Picture1.PSet (-0.1, -6.5), BackColor
    Picture1.Print "delta G (V)"
ElseIf zoomtog% = 1 Then
    Picture1.PSet (-0.1, -1), BackColor
    Picture1.Print "delta G (V)"
End If
'y-axis tick marks
If zoomtog% = 0 Then
    For tmp = -5 To 15
        If tmp Mod 5 Then   'minor tics
            Picture1.Line (-3, tmp)-(-2.95, tmp)
            Picture1.Line (2, tmp)-(1.95, tmp)
        Else
            'major tics
            'right axis
            Picture1.Line (-3, tmp)-(-2.875, tmp)
            Picture1.PSet (-3.1, tmp), BackColor
            out$ = Format(10 ^ tmp, "#.0 E-##")
            Picture1.Print out$
            'left axis
            Picture1.Line (2, tmp)-(1.875, tmp)
            Picture1.PSet (2.3, tmp), BackColor
            Picture1.Print tmp
        End If
    Next
ElseIf zoomtog% = 1 Then
    For tmp = 0 To 15
        If tmp Mod 5 Then   'minor tics
            Picture1.Line (-2, tmp)-(-1.95, tmp)
            Picture1.Line (1, tmp)-(0.95, tmp)
        Else
            'major tics
            'right axis
            Picture1.Line (-2, tmp)-(-1.875, tmp)
            Picture1.PSet (-2.1, tmp), BackColor
            out$ = Format(10 ^ tmp, "#.0 E-##")
            Picture1.Print out$
            'left axis
            Picture1.Line (1, tmp)-(0.875, tmp)
            Picture1.PSet (1.3, tmp), BackColor
            Picture1.Print tmp
        End If
    Next
End If
'Print y-axis labels
If zoomtog% = 0 Then
    Picture1.PSet (-3.5, 7), BackColor
    Picture1.Print "k (/s)"
    Picture1.PSet (2.85, 7), BackColor
    Picture1.Print "log k"
ElseIf zoomtog% = 1 Then
    Picture1.PSet (-2.5, 7), BackColor
    Picture1.Print "k (/s)"
    Picture1.PSet (1.85, 7), BackColor
    Picture1.Print "log k"
End If

tg% = 0
ddG = 0
lastdG = 0
lastval1 = 0
lastval2 = 0
overslope = 0

'Calculate Occupancy of semquinone
K1 = 10 ^ (deltaE / RTmV)
Occup = K1 / (1 + K1)

'Marcus curves

'Calculate and plot a curve for the first electron.

'We calculate the rate constant at different values using an activation energy.
'For Moser-Dutton (MD) curves, the activation energy is that for electron transfer,
'using either the Moser-Dutton algorithm (factor 3.1 contains quantum-mecchanical factor)
'or the classical Marcus term (factor is 4.228 at 298 K), both with deltaG from -zFdeltaE,
'For Bronsted offset curves (MDB or MB), the activation barrier is the sum of deltaG(proton),
'given by the Bronsted term, and deltaG(electron) given using either the Moser-Dutton or Marcus algorithm,
'but with the fraction of deltaG not provided by the Bronsted term
'The curve plotted is for the electron transfer component, and
'this is offset on the logk axis from the standard curve by the Bronsted contribution.

For deltaG = -3.5 * lamda1 To 1.5 * lamda1 Step 0.01  '10 mV steps, deltaG in V
    If MD = True Then                               'If using Moser-Dutton curve
        logk1 = 15 - 0.6 * dist1 - (3.1 * 298 / TK) * ((deltaG + lamda1) ^ 2 / lamda1)
    ElseIf MDB = True Then                          'If using Moser-Dutton curve with Bronsted offset
        logk1 = 15 - 0.6 * dist1 - (3.1 * 298 / TK) * ((deltaG + lamda1) ^ 2 / lamda1) - (pKQH2 - pKISPox)
    ElseIf Marcus = True Then                       'If using Marcus curve
        logk1 = 13 - (beta * (dist1 - 3.6)) / Log(10) - (Fac / TK) * (deltaG + lamda1) ^ 2 / lamda1
    ElseIf MB = True Then                           'If using Marcus curve with Bronsted offset
        logk1 = 13 - (beta * (dist1 - 3.6)) / Log(10) - ((Fac / TK) * (deltaG + lamda1) ^ 2 / lamda1) - (pKQH2 - pKISPox)
    End If
    
            If tg% = 0 Then
            If logk1 >= -5 Then
                If deltaG >= -3 Then
                    Picture1.PSet (deltaG, logk1)
                    tg% = 1
                    lastval1 = logk1
                    lastdG = deltaG
                End If
            End If
        ElseIf logk1 >= -5 Then
            If deltaG >= -3 Then
                EtpRate = 10 ^ logk1
                    If EtpRate > 1.5 * 10 ^ 3 Then
                        LineCol = QBColor(12)
                    Else
                        LineCol = ForeColor
                    End If
                Picture1.Line -(deltaG, logk1), LineCol
                'Find change in logk for this 10 mV step
                dlogk1 = logk1 - lastval1
                ddG = deltaG - lastdG
'We print the rate and slope on the graph,
'using the intercept of the Marcus curve with the value for dE,
'or for dEelect, depending on whether we are using a Bronsted term.
                If MD = True Or Marcus = True Then
                    If (-deltaE1 > lastdG And -deltaE1 < deltaG) Or -deltaE1 = deltaG Then
                        slope_logk1 = dlogk1 / 10
                        Picture1.PSet (deltaG, logk1)
                        Picture1.ForeColor = QBColor(12)
                        Picture1.Print Format(10 ^ logk1, "#.##e-##"); "/s, "; Format(slope_logk1, "##.####") & "/mV"
                        Picture1.PSet (deltaG, logk1)
                        Picture1.ForeColor = QBColor(0)
                        Rate1st = 10 ^ logk1
                    End If
                ElseIf MDB = True Or MB = True Then
                   If (-deltaE1elec > lastdG And -deltaE1elec < deltaG) Or -deltaE1elec = deltaG Then
                        slope_logk1 = dlogk1 / 10
                        Picture1.PSet (deltaG, logk1)
                        Picture1.ForeColor = QBColor(12)
                        Picture1.Print Format(10 ^ logk1, "#.##e-##"); "/s, "; Format(slope_logk1, "##.####") & "/mV"
                        Picture1.PSet (deltaG, logk1)
                        Picture1.ForeColor = QBColor(0)
                        Rate1st = 10 ^ logk1
                    End If
                End If
                
                lastval1 = logk1
                lastdG = deltaG
            End If
        End If
Next
tg% = 0
overslope = 0
'Calculate and plot a curve for the second electron
For deltaG = -3.5 * lamda2 To 1.5 * lamda2 Step 0.01  '10 mV steps
    
    If MD = True Then
        logk2 = 15 - 0.6 * dist2 - (3.1 * 298 / TK) * (deltaG + lamda2) ^ 2 / lamda2
    ElseIf MDB = True Then
        'deltaGe = deltaG - (pKQH2 - pKISPox)
        logk2 = 15 - 0.6 * dist2 - (3.1 * 298 / TK) * (deltaG + lamda2) ^ 2 / lamda2
    ElseIf Marcus = True Then
        logk2 = 13 - (beta * (dist2 - 3.6)) / Log(10) - (Fac / TK) * (deltaG + lamda2) ^ 2 / lamda2
    ElseIf MB = True Then
        'deltaGe = deltaG - (pKQH2 - pKISPox)
        logk2 = 13 - (beta * (dist2 - 3.6)) / Log(10) - (Fac / TK) * (deltaG + lamda2) ^ 2 / lamda2
    End If
   If tg% = 0 Then
            If logk2 >= -5 Then
                If deltaG >= -3 Then
                    Picture1.PSet (deltaG, logk2)
                    tg% = 1
                    lastval2 = logk2
                    lastdG = deltaG
                End If
            End If
        ElseIf logk2 >= -5 Then
            If deltaG >= -3 Then
                EtpRate = Occup * 10 ^ logk2
                    If EtpRate > 1.5 * 10 ^ 3 Then
                        LineCol = QBColor(14)
                    Else
                        LineCol = ForeColor
                    End If
                Picture1.Line -(deltaG, logk2), LineCol
                'Find change in logk for this 10 mV step
                dlogk2 = logk2 - lastval2
                If (-deltaE2 > lastdG And -deltaE2 < deltaG) Or -deltaE2 = deltaG Then
                    slope_logk2 = dlogk2 / 10
                    Picture1.PSet (deltaG, logk2)
                    Picture1.ForeColor = QBColor(14)
                    Picture1.Print Format(Occup * 10 ^ logk2, "#.##E-##") & "/s, "; Format(slope_logk2, "#.####") & "/mV"
                    Picture1.PSet (deltaG, logk2)
                    Picture1.ForeColor = QBColor(0)
                    Rate2nd = Occup * 10 ^ logk2
                End If
                'savelast = lastval1
                lastval2 = logk2
                lastdG = deltaG
lastval2 = logk2
                lastdG = deltaG
            End If
        End If
        If overslope = 0 Then
            drv = Abs(dlogk2 / 10)       'dlogK/mv
            If drv >= 0.0095 Then
                Picture1.PSet (deltaG, logk2)
                Picture1.ForeColor = QBColor(14)
                'Picture1.Print Format(dlogk2 / 10, "#.####")
                Picture1.PSet (deltaG, logk2)
                Picture1.ForeColor = QBColor(0)
                overslope = 1
            End If
        End If
Next

If MD = True Or Marcus = True Then
    Picture1.Line (-deltaE1, -3)-(-deltaE1, 12.5)
    Picture1.PSet (-deltaE1, 13), BackColor
    Picture1.Print "E1"
ElseIf MDB = True Or MB = True Then
    Picture1.Line (-deltaE1elec, -3)-(-deltaE1elec, 12.5)
    Picture1.PSet (-deltaE1elec, 13), BackColor
    Picture1.Print "E1'"
End If
Picture1.Line (-deltaE2, -3)-(-deltaE2, 13)
Picture1.PSet (-deltaE2, 13.5), BackColor
Picture1.Print "E2"
Picture1.Line (deltaG_bound, -3)-(deltaG_bound, 13.5)
Picture1.PSet (deltaG_bound, 14), BackColor
Picture1.Print "Gb"
Picture1.Line (deltaG_overall, -3)-(deltaG_overall, 14)
Picture1.PSet (deltaG_overall, 14.5), BackColor
Picture1.Print "Go"
Picture1.DrawStyle = 1
Picture1.Line (-3, Log10(1500))-(2, Log10(1500))
Picture1.DrawStyle = 0
End Sub


Private Sub CalcLamda1Cutoff()
oldGact = 0
For lamda = 0.2 To 4 Step 0.001
    Gact = ((-deltaE1 + lamda) ^ 2) / (4 * lamda)
    If Gact >= testEact1 And oldGact < testEact1 Then
        Lamda1cutoff = lamda
        Text8.Text = Lamda1cutoff
    End If
    oldGact = Gact
Next
If lamda1 > Lamda1cutoff Then
    Text1.BackColor = QBColor(11)
Else
    Text1.BackColor = QBColor(15)
End If
End Sub

Private Sub CalcLamda2Cutoff()
oldGact = 0
For lamda = 0.1 To 4 Step 0.001
    Gact = ((-deltaE2 + lamda) ^ 2) / (4 * lamda)
    If Gact >= testEact2 And oldGact < testEact2 Then
        Lamda2cutoff = lamda
        Text9.Text = Lamda2cutoff
    End If
    oldGact = Gact
Next
If lamda2 > Lamda2cutoff Then
    Text6.BackColor = QBColor(11)
Else
    Text6.BackColor = QBColor(15)
End If
End Sub


Private Sub CalcGvals()

'Calculate deltaG values for partial reactions

'Overall, where deltaG values are in V (deltaG / F)
'From deltaG/F = -sigma(zE), with E donor couple negative.
deltaG_overall = -((EmISP + EmbL) - 2 * EmQQH2) / 1000
deltaG_bound = -((EmISP + EmbL) - 2 * EmQQH2bound) / 1000

'Set value for deltaE for 1st electron transfer
deltaE = HScroll5.Value 'value in mV
deltaE1 = deltaE / 1000  'value in V

'Calculate fraction due to H+ transfer
deltaGprot = RTV * (pKQH2 - pKISPox)      'Bronstedt term, in volts

'Fraction due to electron transfer
deltaE1elec = deltaE1 + deltaGprot          'plus, because dG = -zfdE

'Em QHdotQH2 (in mV)
EmQHdotQH2 = EmISP - deltaE

'Em QQHdot (in mV)
EmQQHdot = 2 * EmQQH2bound - EmQHdotQH2

'Ks
logKs = (EmQQHdot - EmQHdotQH2) / RTmV
Text10.Text = Format(logKs, "##.##")

'2nd electron transfer
deltaE2mV = (EmbL - EmQQHdot)       ' (in mV)
deltaE2 = (EmbL - EmQQHdot) / 1000
'MsgBox "E1 = " & Str$(deltaE1) & ", E2 = " & Str$(deltaE2) & ", Gb = " & Str$(deltaG_bound * 1000)
If deltaE + deltaE2mV <> -deltaG_bound * 1000 Then
    MsgBox "Error" & "E1 = " & Str$(deltaE1) & ", E2 = " & Str$(deltaE2) & ", Gb = " & Str$(deltaG_bound * 1000)
End If

End Sub

Private Sub Text10_Click()
Fix_Ks% = 1
Fix_deltaE% = 0
Text10.BackColor = QBColor(10)
Text5.BackColor = QBColor(15)
End Sub

Private Sub Text5_Click()
Fix_Ks% = 0
Fix_deltaE% = 1
Text10.BackColor = QBColor(15)
Text5.BackColor = QBColor(10)

End Sub