- Email: [email protected]
Electrostatic potentials in membrane systems
Bulletin of Mathematical Biology, Vol. 45, No. 1, pp. 69-90, 1983. Printed in Great Britain.
0092-8240183]010069-22503.0010 Pergamon Press Ltd. © 1983 Society for Mathematical Biology
E L E C T R O S T A T I C P O T E N T I A L S IN MEMBRANE SYSTEMS • JACEKT. DUNIECand SIDNEYW. THORNE CSIRO Division of Plant Industry, P.O. Box 1600, Canberra 2601, Australia
A membrane with an arbitrary distribution of fixed charges inside and on its surfaces is considered. A procedure for calculating the local electrostatic potential at an arbitrary point of the system is described and its validity discussed. This procedure is based on the linearization of the 3-dimensional Poisson-Boltzmann equation around an exact 1dimensional solution.
1. Introduction. Several phenomena at a molecular level in biological membranes appear to be influenced by local electric fields. This is evident in particular in electron transfer reactions; for example, the redox reactions mediated by the cytochrome c depend significantly on the distribution of charges on this protein (Speck et al., 1979). Also, in photosynthetic membranes local transmembrane electric fields are likely to provide the necessary bias for the vectorial charge separation induced by light (Duniec and Thorne, 1979, 1980, 1981). Other examples are discussed in the review by Waltz (1979). It seems important, therefore, to be able to predict the magnitude of such fields in biomembranes as a function of charges inside and on the surface of the membrane as well as a function of ionic conditions in the surrounding solution. The simplest approach is based on the GouyChapman theory of electric potentials in electrolytes near charged surfaces of colloidal particles (Barlow, 1970). In this model one assumes that the membrane can be considered as an infinite flat layer sandwiched between the two half-spaces occupied by aqueous solutions. The surfaces of the membrane are charged with uniform surface densities which are representative of real densities of (discrete) fixed charges. This approach is adequate for the description of phenomena which depend on average electric potentials and fields. Examples can be: calculations of electrostatic repulsive force between the membranes (Duniec et al., 1979), local ionic concentrations near the membrane surface (Barber et al., 1977) and electrophoretic mobility of membranous particles (Nak69
70
J . T . DUNIEC AND S. W. THORNE
atani et al., 1978). However, for the processes at the molecular level the local electric fields and potentials are important, thus a model which accounts for the discrete nature of fixed charges and particular dielectric properties of the water-membrane interface is needed. Several approaches to build the model for the calculations of the discrete charge effects on the electric potentials in membranes have been made in the past [Sauv~ and Ohki (1979) and references therein]. The most recent method, developed by Sauv~ and Ohki (1979), allows for the calculation of electric potentials at the water-membrane interface for an arbitrary distribution of fixed charges on the membrane surface. Their procedure is based on the linearization of the Poisson-Boltzmann equation (describing electric potentials in electrolytes) around a certain 'smeared surface charge' solution. The Hankel transform technique is then applied to the resulting system of linear equations. A similar method was previously used in problems of adsorption of ions at the metal-water interface by Loeb (1951) and Levine et al. (1972). The work of Sauv6 and Ohki (1979) deals with a single water-membrane interface while the membrane is considered infinitely thick. Such an approximation is valid for calculations of local potentials due to charges near the surface of the membrane as there is only weak coupling across the membrane (Duniec and Thorne, 1981). The purpose of this paper is to provide a mathematical method for calculating local potentials anywhere in the membrane system. The fixed charges may be positioned on each side of the membrane as well as inside it. Such a method is needed for the problems related to the electron transfer reactions inside and across the membrane, where the internal charges and finite membrane thickness must be considered. An example of application of this method has already been published (Duniec and Thorne, 1981).
2. Mathematical Preliminaries. The fundamental physical law underlying the electrical phenomena in membrane systems is the Gauss' law of electrostatics: V. (e(r)V~(r)) = - 47rp(r),
(1)
where ((r), ~(r) and p(r) are the dielectric permittivity, electric potential and external charge density, respectively, at a point r = (x, y, z). The operator V = ((O/ax),(O/ay),(a/az)) and (.) denotes the scalar product. Throughout this work we shall be concerned with the system consisting of (see Figure 1): (1) th e aqueous region extending from z = - o o to z=-(a+d) andfrom z=a+dto z = + ~ ; ( 2 ) two polar regions: - a -
ELECTROSTATIC
POTENTIALS
ol° xy
I I I I
Wa~er
IN MEMBRANE
o o
SYSTEMS
71
I I
o
I o% I
o °
Membrane o
a 7-
Io
o 2d~
I o I I
----~ o
%
I
L I
o o
o
o
I I
o
I I
3
O
o
,:9
o
O
I
Figure 1. The coordinate system is chosen with the origin in the middle of the membrane and X, Y axes in its plane. The water region surrounds the membrane from both sides. Its interfaces with the membrane define the polar regions of thickness a. The dots in the membrane and the polar regions denote discrete charges.
d < z < -7 d and d < z < a + d; and (3) the hydrophobic portion of the membrane - d < z < d. The dielectric permittivities of each region will be assumed to be constant. Let us denote the dielectric constants as follows: ew for the water region, ep for the polar region and Em for the hydrophobic portion of the membrane. Let us first consider the aqueous region. We will restrict our attention to the case w h e n there are only salts of mono- and divalent cations dissolved in the water. All ionic species contribute to the volume charge density which, according to the Boltzmann's law, depend on the electric potential: p(r) = q[ - C - exp (qW/kT) + C + exp ( - q ~ / k T ) + 2C 2+ exp ( - 2q~/kT)],
(2) where q is the elementary charge, and C ÷, C - and C 2+ are the bulk concentrations of m o n o v a l e n t cations and anions and divalent cations respectively. F o l l o w i n g the notation of Ninham and Parsegian (1971), w e introduce the dimensionless electric potential tO= q'q/kT and the parameters 0 ~< rt ~< 1 and no such that: C 2+ = (1/2)~qno and C- = no. Since
72
J . T . D U N I E C A N D S. W. T H O R N E
the bulk solution must be neutral, i.e. 2 C 2+ + C + - C = 0, we can see that C + = (1 - "0)no. N o w f r o m (1) and (2) we obtain the following equation for tO(r): VztO = (/(2/2)[exp (tO) - (1 - "0) exp ( - tO) - "0 exp ( - 2tO)], or a + d < z < + o ~ .
-oo
(3)
The Debye's c o n s t a n t / ( is defined by /(2= 8.rrnoq2/ewkT. Equation (3) is a semilinear elliptic partial differential e q u a t i o n of the second order to which general t h e o r e m s of existence and u n i q u e n e s s of solutions to boundary value p r o b l e m s can be applied (Gilbarg and Trudinger, 1977). Practical m e t h o d s for solving such equations are, h o w e v e r , complicated. In A p p e n d i x 1 we shall discuss sufficient conditions for which the procedure based on the linearization of (3) a r o u n d an a priori known solution of the same equation is valid. Let tOo(Z) be a solution to (3), i n d e p e n d e n t of x and y, such that it vanishes at infinity. L e t
(4)
tO,,w(r) = tO(r) - tOo(Z).
Substituting (4) into (3), expanding (3) in terms of tO~,wand retaining only the leading terms, we obtain
2
2 '~.
v tO,, =/(o(Z).tO,.w,
(5)
where [in notation of Sauv6 and Ohki (1979)] Ko(Z) = (K2/2)[exp (tOo(Z) + (1 - n) exp ( - tOo(Z)) + 2"0 exp ( - 2tOo(Z))]
(z<-a-d
or z > a + d ) .
(6)
For notational c o n v e n i e n c e let us introduce the characteristic f u n c t i o n s for the polar regions and the interior of the m e m b r a n e :
X~(Z)
=fl'
t
X,,(z) X+ p(z)=
for - d - a < z < - d 0, elsewhere
=,[1, for d < z < d t 0, elsewhere
( 1,
ford
..
(7a)
(7b)
(7c)
E L E C T R O S T A T I C P O T E N T I A L S IN M E M B R A N E SYSTEMS
73
In addition, we shall use the superscript ( - ) for all the functions defined for z < d and the superscript ( + ) for those defined for z > d. As already implied in equations (4) and (7), the subscripts p, m and w will denote the quantities referring to the polar regions, membrane and water respectively. Let us assume that there are a total of N fixed charges at positions R~, and of valences v~ located in the membrane and the polar regions. Using equation (1) and the introduced notation, we can write the Poisson equation for the dimensionless potentials to~--in the polar regions and tom--in the membrane: N
V2to; = - ( 4 ¢ r q 2 / e p k T ) ~
v~8(r - r~)x~(z~)
(8)
v i S ( r - ri)xm(zj.
(9)
i=l
N
V2tom = - ( 4 7 r q 2 [ e m k T ) ~ i=l
The boundary conditions that are imposed on the potentials O?w, to~ and 0m can be derived from: (1) the continuity of the total electrostatic potential to(r) and the quantity e(dto/dz) (where E may be ew, Ep or Era) across the boundaries, and (2) the electroneutrality of the bulk solution which is equivalent to the requirement that the electrostatic potential vanished for z ~ +oo. These boundary conditions can be expressed as follows: lim O~,w(r)= 0
(lOa) +
too( +-- a +- d ) + to~,w(O, ++-a +- d ) = to;(O, +- a +- d )
Ew
+ Ew z=+-a±d
C
= Ep z=±
a±d
to;(p, _+ d) = to (p, -+ d)
(lOb)
z=±
a±d
(lOd) (lOe)
where p = (x, y). Linear equations (5, 8, 9) with the boundary conditions (10) can be solved by representing the potential in each region as the sum of contributions from individual charges (Sauv6 and Ohki, 1979) plus an additional 1-dimensional component due to a possible asymmetry between the potentials too(Z) and t0o(Z). Thus we seek the solution to equations (5, 8, 9) in the form:
74
J. T. DUNIEC AND S. W. THORNE N
O~,w(r) = ~bw(Z) + ~, ~w,~(r; r,)
(lla)
,=1
N
qJ;(r) = 6~(z) + Y. ~.,(r; r,)
(llb)
i=1
N
tOm(r) = q~,,(z) + ~ ~m,~(r; r,),
(llc)
i=l
where ~;,~, g~,~ and ~,~ are the individual charges contributions to the potentials 4'~.w, q'; and 4'm respectively, and qQ, 4~; and 4~,, are introduced to match the inhomogeneous boundary conditions (10b, c). The potentials ~,,, ~,~ and ~,~ obey the equations (i = 1, 2 , . . . , N): 2
m
2
±
Vr~,~(r; ri) = Ko(Z)~w,i(r; r,) 2
--+
--
4"a'q 2
(12a) +
V,~:p,,(r; r,) - - e p k T v~8(r - ri)x~(z~) Vr2G~.,(r; r,) -
4~rq 2 7~v,8(r-
r,)xm(z3.
(12b) (12c)
These potentials satisfy the following boundary conditions: ~,i(p, -+ a --- d; ri) = ~,~(p, -+ a -+ d; r~)
0G,,_ aG,,
ew0---7- - e p - ~ - , z = -+ a -+ d ~.,(p, ___d; ri) = ~m,,(p, -+d; ri)
(13a) (13b) (13c)
eo O~;'t aG,, oz = e m ~ , z = +- d
(13d)
lira ~w,~(r; ri) = 0.
(13e)
On the other hand, the potentials rbw(Z), 49;(z) and ~b,,(z) appearing on the r.h.s, of equation (11) obey the equations: ~ =
~:~(z)q~;
(14a)
ELECTROSTATIC POTENTIALS IN MEMBRANE SYSTEMS
d%; _d%m - ~ =0
75
(14b)
with the boundary conditions which depend on the choice of too(Z): &~( --- a + d) + too( ± a ± d) = ~b~(± a --- d) ew
+ew
=ep
, z=+a+d
&~(--- d) = 6•(± d) ep ~ z
(15a)
(15c)
~ =era dd4~ z' , z = ± d
(15d)
lira & ; ( z ) = 0.
(15e)
We have thus divided the linear system of equations (5, 8, 9) with the boundary conditions (10) into N independent systems (12, 13), each corresponding to an individual fixed charge in the polar region or inside the membrane, and one additional system (14, 15) whose solutions must be added to the sum of contributions of individual charges so that the total potential satisfies boundary conditions (10). Note that although boundary conditions (13) do not involve the function too(Z), the contributions of individual charges in each region depend on the choice of t0o(Z) through the function K~(z) [cf. eq (6)1 appearing in (12a). 3. C o n s t r u c t i o n o f t h e S o l u t i o n s to t h e L i n e a r i z e d P r o b l e m . We proceed now to construct the solutions to equations (12, 13) and (14, 15). Note that by symmetry, the potentials of individual charges depend on the absolute value of the difference p - p~, so that we can write:
•
~L(o, z;oi, z,) = ~;,,(R, z; z,)
(16a)
C A p , z; p,, z,) = C A R , z; z,)
(16b)
~,.,~(p, z; Oi, z,) = ~,.,~(R, z; zi),
(16c)
where
R=lp-p,I.
76
J.T.
DUNIEC
A N D S. W . T H O R N E
Accordingly, we can rewrite the operator V~ as: d 2
2
0 2
V~z = V~ + d-~, where VR = ~-~ + R -~ O0R '
(17)
and substitute equations (16, 17) into (12). In the case of (12a) this yields: /
I
J ~2x
n,~raz
I
(18)
= Ko(Z)~;.,(R, z ; z,).
"
Let us apply the Hankel transform defined as: dr
(19a)
f ( R ) = I o a f ( a ) S o ( a R ) da
(19b)
= Io Rf(R)Jo(aR)
to equation (18). After some algebra we arrive at: d. 2 ~+ 2 ~± d z 2 ¢~.i(a, z; zi) = (ct 2 + Ko(Z))~w,i(ol , z ; Zi).
(20)
Before applying a similar treatment to (12b, c) let us represent the functions ~.i and ~m.i as: 2 +
+
viq
+
~p.~(R, z; zi) = ~,,(R, z, z~) + epkTv,(R~ + (z - z~):) Xp(z~)
(21a)
2
+- z ')" ~m,i(R, z; z~) = ~,~(R, z, z~) + ~ , . k T v , ( R v~q ~ + (z - z~)~) X'n(
(21b)
Substituting (21) to (12b, c) and bearing in mind (16 and 17), we arrive at the following equations for ~.~ and ~m.~:
(
2
V2+
~p.~(R,z; z~) = 0
az /
"
Now the application of the Hankel transform to (22) yields:
(22a)
ELECTROSTATIC POTENTIALS IN MEMBRANE SYSTEMS
77
2~+
d ~19,i 2 "'+ dz 2 = a ~,~
(23a)
d2~,~,i 2dz = c~ ~,,,,i.
(23b)
Equations (20, 23) can be solved in a general form. In particular, the solution to (20) which vanishes at infinity can be written as [using the notation of Sauv6 and Ohki (1979)]:
~w,, = A~,iQ+-(z, a, "q) exp ( +- Paz),
(24)
where A~,, are certain constants to be found from the boundary conditions (13), and P is defined below: P = (1 + K2(1 + rl)/a2) '/2.
(25)
The functions Q-+(z, a, r/) (Sauv6 and Ohki, 1979) are described in Appendix 2. The general solutions to (23) can be written as:
~,, = A~,, exp ( - az) + B ~,, exp (az)
(26a)
~m,, = A,~,i exp ( - az) + B,,,., exp (az).
(26b)
In order to determine the 8 constants A~,,; A~.i; B~,i; Am,i; and Bin,i, we will derive a system of algebraic equations from the boundary conditions (13a-d). To this end, we use representations (16) and (21) in (13a-d) and then apply the Hankel transform (19) with respect to R = ]p -Oil. Using the well-known identity:
f R V/(R2 + 1(z -
zi) 2)
Jo(aR) d R
= 1 exp ( -
alz - z,I),
(27)
the above-described procedure yields: ~ , , ( a , ---a - d) = gp,i(a, -- a +- d) + /3, exp ( - ala + d -T-zil)X~(zi) ~ep (28a) d~w,, =
Ew d z
Ep
dd~/3,
exp ( - ala + d w z,l)Xp(z, )
(28b)
78
J . T . D U N I E C A N D S. W. T H O R N E
atz=
+_a+-d,
~;.~(a, +- d) + X~(Z,)a-~ exp ( - a]d -T z,]) = ~.,.,(a, +_ d) + X.,(z~) a ~ ~ exp ( - a ] d w- z~[) (28c)
e p ~ - - - / 3 i exp ( - aid -T- zil)X~(zi) = d~m.~ _
e,,--d~- + [3i exp ( - aid -~ zil)X,.(zi) (28d) at z = -+ d, where [3i = v~q2/kT. Note that equations (28) f o r m a s y s t e m of eight equations since all of the above formulae can be read with the u p p e r sign and with the lower sign. Let us introduce the eight-dimensional vector X consisting of the coefficients appearing in equations (24, 26): +
+
+
(X') r =(A~,,;A~,,;B~,,;Am,~;Bm,,;A~,,;Bp,~;Aw,~).
(29)
Let Y~ be the 8-dimensional vector w h o s e c o m p o n e n t s Y~, ] = l, 2 . . . . . 8 are defined as follows:
YI =
(13,1e,,a)x;(z,)
exp ( -
Y~2 = fl,X;(Z,) exp ( -
ala +
ala + d + z,I)
(30a)
d + z,I)
(30b)
Y~ = (fl,/a)(Xm(z,)/e,,, - X;(z,)/e~,) exp ( - aid + z,I)
(3Oc)
Y4 = (,Sffa)(X.,(z~) + X-p(zO) exp ( - aid + z,[)
(30d)
Y5 = ([3,/a)(xS(z,)/ep
- X . . ( z , ) / e . . ) e x p ( - aid - z,I)
(30e)
-b
Y6
=
(fl,la)(x,,(z,) + x,,,(z,))
exp
( -
a i d - z,I)
(30f)
Y~, = -(fl,/e~a)x;(z,) exp ( - ala + d - z,I)
(30g)
Y ; = ~,x+~(z,) exp ( - a l a + d - z,I).
(30h)
With the above defined notation, substitution of equations (24, 26) into
ELECTROSTATIC POTENTIALS IN MEMBRANE SYSTEMS
79
equation (28) results in a system of algebraic equations for the coefficients XI(j = 1 , 2 , . . , 8 ) which can be written in the following concise form: A X i = yi,
(31)
where the matrix A is defined by: p q 0 0 0 0 0 0
-
- -
-s
S
--1
r
A=
aeps t - ept 0 0 0 0
- a~ps t -1 ept -~ 0
-1
0
0
0
0
0
0
0
0
0
0
0
0
0
-
t
-
t -1
e,~t t -~
- e,~t -~ t
0
-- e,nt -1
emt
0 0
0 0
0 0
0 - t -~ ept ~ s -~ - epas -~
0 -t -- ept s %cts
0 0 0 p+ q+
(32) with the following notation: p-*(a) = 7 Q-*( -- a -- d, a, 7/) exp ( - P a ( a + d))
[ OQ*- 11
(33a)
-7- pQ*-( +_ a + d, a, ")In exp ( - P a ( a + d))
(33b) s ( a ) = exp ( a ( a + d))
(33c)
t ( a ) = exp ( a d ) .
(33d)
Upon solving (31) we obtain the coefficient vector Xi(a) [cf. (29)]. Thus we can find the Hankel transform of the potentials ~ , i ( a , z ; z~), ~,i(a, z; z~) and g,,,i(a, z; z~) from equations (24, 26). Subsequently, by inverting the Hankel transforms according to formula (19b) and using (21), we can find the ith charge contribution to the total electric potential in each region. In actual numerical calculations one needs to deal with the fact that for a = 0 the determinant of the matrix A [equation (32)] vanishes. However, the transforms ~;,~(a,z;zi), ~ . i ( a , z ; z O and ~,,,~(a, z; zi) are well defined for a = 0, which can be checked by examining the following limits: ~;,~(0, z; z~) = lim (A;,~(a)Q*-(z, a, "q) exp ( ~ P a z )
(34a)
80
J. T. DUNIEC AND S. W. THORNE
1
(;,i(~,z; zi) = lim Az,i(a)exp ( - az) + B z,i(a)exp (az) a-0
21.
+-J-exp
( - a i r - zil)x;(zi)]
E P ~
(34b)
{m,i(0, Z ;zi) = lim Am,i(a)exp ( - az) + Bm,i(a)exp (az) u -0
where A:,i, A;,i, B;,i, Am,i and BmViare the components of the vector x i ( a ) [equation (29)l. Performing the necessary algebra and using the notation
we obtain:
x
Q-(z7 O, ')
exp ( ~ ' ( 2+ a + d))
Q-( - a - d, 0, TJ)
(35c)
E L E C T R O S T A T I C P O T E N T I A L S IN M E M B R A N E S Y S T E M S
-aT~,-+~,,(0, ~l
81
z; z,) = (p -q ?P+_P +q-)(q-O + p-)
Q+(z, o, n)
x Q+(a + d, 0, a~) exp ( - K'(z - (a + d))).
(35g)
For the construction of the full solution to equations (8-10) we n e e d to find the potentials ~bw(Z), ~b~(z) and ~bm(z) appearing on the r.h.s, of equations (11) and satisfying equation (14) with the b o u n d a r y conditions (15). In general, these functions can be r e p r e s e n t e d as follows: ¢~,(z) = A~Q+-(z, 0, ~) exp ( ~ K'z) +
(36a)
&~(z) = A ~ z + B y
(36b)
(~m(Z) = A m z + Bin.
(36c)
Introducing the v e c t o r +
+
X ° = ( A : , A~, By, Am, Bin, Ap, Bp, Aw),
(37)
and using notation (33a, b), we can write b o u n d a r y conditions (15) as follows: --
O
0
0
p X 1 + X 2 ( a + d) - X 3 = - a -
0
0
--
(38a)
q X1 - epX2 = - b -
(38b)
- d X ° + X ° + dX°4 - X ° = 0
(38c)
epX ° - E~X°4 = 0
(38d)
dxo4+ X 5o - d X 6o -
(38e)
X ° = 0
e m X ° - epX ° = 0 +
0
p X8+(a+d)X°6+X +
0
0 6 =
(38f) °=a
+
+
(38g)
b ,
(38h)
a ± = tOo( -+ (a + d))
(39a)
q
X8
"~- E p X
w h e r e use was m a d e of the notation:
82
J . T . DUNIEC AND S. W. THORNE
b -+ =
ew
z = - (a + d).
(39b)
Once the solutions to (38) have been found, the full solution to the linearized problem can be obtained. As we have shown, this solution is determined by (a) the function too(Z), which is an exact solution to the 1-dimensional nonlinear Poisson-Boltzmann equation, and (b) the positions and valences of individual charges within the polar regions and the membrane. Our approach of approximating the solution to the 3-dimensional nonlinear equation (3) by a sum of an exact 1-dimensional solution to (3) and a solution to 3-dimensional linear equation (12a) is a first step in the Newton-Kantorovich method of constructing solutions to nonlinear problems (Krasnosel'ski e t a l . , 1969). In Appendix 1 we find the sufficient conditions under which this procedure as applied to the nonlinear Poisson-Boltzmann equation is convergent. Now we are concerned with finding an appropriate 'zeroth order' approximation t0o(z). An obvious choice is to let too(Z) be the solution to the 1-dimensional nonlinear Poisson-Boltzmann equation with the following boundary condition at z = - (a + d) (Sauv6 and Ohki, 1979): b ~- = ew
-
S
vi,
(40)
_+
where S is the area of the membrane, and I_+ are sets of indices denoting the charges in the r.h.s. ( + ) and the 1.h.s. ( - ) polar regions. Although this choice is not optimal, as there is some degree of coupling across the membrane, the possible discontinuity of the electric potential is rectified in the next step, i.e. the first-order approximation. We conclude this section by outlining a computational procedure which can be used in practical calculations. 1. Solve 1-dimensional Poisson-Boltzmann equation (3) with the boundary conditions (40). Obtain too(z). 2. Solve system of equations (38) using too(Z). Obtain ~bw(Z), $~(z) and qSm(z) according to (36, 37). 3. For each a needed to determine numerically the Hankel transform, and for each i, find the coefficient vector X i (29) from (31). Use an appropriate subroutine for solving a system of linear equations rather than explicit forms for X i, since these are very complex. 4. Using the coefficient vector X ~ (29) find the potentials ~,~, ~,i and ~,,.i by applying the inverse Hankel transform to equations (24, 26). For a ~ 1, use the asymptotic forms expressed by (35).
ELECTROSTATIC POTENTIALSIN MEMBRANESYSTEMS
83
5. The value of the total potential at a particular point r = (p, z) can be found by adding appropriate solutions obtained in steps 1-4. For example, the potential at the 1.h.s. polar region will be: N
4 (o, z) = 6 ; ( z ) + 5". ;(Io -o 1, z; z,) i=1
etc. 4. C o n c l u d i n g R e m a r k s . The purpose of this paper is to present a c o m p r e h e n s i v e description of the m e t h o d of calculating the local electric potentials in m e m b r a n e systems. We have avoided giving a specific example on purpose, since one application of the m e t h o d has already b e e n published (Duniec and Thorne, 1981). In this r e f e r e n c e the effects of salt conditions, positions of charges and dielectric properties of a w a t e r m e m b r a n e interface on the local potentials w e r e discussed. N e v e r t h l e s s it is worthwhile to m e n t i o n here the possible usefulness of the m e t h o d in our u n d e r s t a n d i n g of electron transfer processes in m e m b r a n e s . Electron transfer rates are affected by the electrostatic potential difference b e t w e e n the points w h e r e the donor and the acceptor reside. A change in this quantity m a y accelerate or slow d o w n the reactions. As we have shown r e c e n t l y (Duniec and Thorne, 1981), changes of the value of the dielectric permittivity of the polar regions as well as adsorption of ions m a y have a significant effect on the potential difference across the m e m b r a n e . It is possible, therefore, that the control of electron transfer rates in m e m b r a n e s m a y be exercised by such e n v i r o n m e n t a l changes. A n y relevant example would require a detailed k n o w l e d g e of the structure of the m e m b r a n e . As m o r e information on it b e c o m e s available, m a t h e m a t i c a l tools for quantifying electrical p h e n o m e n a at the molecular level will be needed. APPENDIX 1
The linearization of the Poisson-Boltzmann equation applied in this work is the first step in the Newton-Kantorovich method of solving nonlinear problems (Krasnosel'ski et al., 1969). Below we derive sufficient conditions for the applicability of this method to the nonlinear Poisson-Boltzmann equation and thus the validity of the linearization procedure. Consider the equation Fx = 0,
(A1)
where F is a Frechet ditterentiable nonlinear operator defined on a domain Dt of a Banach space B~ and transforming it into a Banach space B2. Let x* be a solution to (A1) and suppose that the nth approximation to x* has been found. Consider now the linear
84
J.T. DUNIEC AND S. W. THORNE
equation:
Fx. + F'(x.)(x - x.) = 0.
(A2)
If the linear operator [F'(x.)] ' exists then the solution to (A2) can be written as
x.+, = x. - [F'(x.)]-'Fx..
(A3)
The above formula defines the iterative Newton-Kantorovich method of solving (A1). The following theorem states sufficient conditions for the convergence of the above iterative process (Krasnosel'ski et al., 1969). THEOREM. Let F be determined in a certain ball S(xo, R) with the centre at Xo and the radius R, let its Frechet derivative F'(x) (treated as an operator defined on S(xo, R) and with values in the Banach space of continuous linear operators transforming B, into BE) satisfy the Lipschitz condition in S(xo, R): LIE'(X)
-
(A4)
f ' ( y ) l l ~ L I I x - YlI,
let Fo = I[f'(xo)-'ll, To = IIF'(Xo)-'FXoll
(AS)
ho = FoLno < 1/2.
(A6)
and
Then the sequence {x,},=o,, .... defined iteratively by (A3) converges to the solution x* E S(xo, R), where R > (1 - (1 - 2ho)'l:)rto/ho.
(A7)
We shall now apply this theorem to the Dirichlet problem of the nonlinear PoissonBoltzmann equation determined in a slab SL ={(x, y, z ) ; 0 < z 0 when the limit d --*~ is taken. Consider, then, the following equation in SL: v % = K2f(g,)
(A8)
with the boundary conditions: 0(rs) = Os(rs), where 0s is a bound, continuous function and r, = (x, y, 0) or (x, y, d). Let ~o(z) be a 1-dimensional solution to (A8), i.e. d20° = ~2f(Oo). dz 2
(A9)
Representing ~ as ~ = 6o + ~1, we can rewrite ( A l l ) as follows: V2~11-- K2(Z)I~I = K2 f(00 ~- ~1)-- f(o0) --
~Jl
--- g(O,),
(AI0)
where we have also defined the function g(t~,). If we let KZo(z)= K2f'(tko(Z)), then
g(o) = g'(O) = 0. Let Go(r, r') be the Green's function for the Dirichlet problem of the equation:
V20, - Ko2(Z)t~,= 0.
(All)
ELECTROSTATIC POTENTIALS IN MEMBRANE SYSTEMS
85
By means of Go(r, r') we can write (A10) in the integral form: tpl(r) = fSL d3r'G°(r' r')g(t~l(r')) - f~=o d2r5 0~Go(r, r's)t~ls(r's) + ~ =a dSrs 0 ~ G ( r , r's)t~,s(r's),
(A12)
where 6~s(rs) = ~ ~s(r~)- 60(0) for z = 0 b k , ( r ~ ) - ~o(d) for z = d" Equation (A12) defines the operator F acting on functions q'l: =- ~sL d3r'G°(r' r')g(t~(r')) + b(r) - t~x = 0.
(13)
The operator F acts in the Banach space of continuous functions C(SL) defined in the slab SL and equipped with the norm:
(A14)
I1~11= sup Iq,(r) I. O
The Frechet derivative of the operator F at some point $~E calculated from (A13):
C(SL)
can be readily
F'(t~l)t~ = fo
(A15a)
It is also apparent that F'(0)~b = - 6
(A15b)
r ° = II(F'(0)-')ll-- 1.
(MSc)
In order to estimate the Lipschitz constant of the operator F we perform the following:
I(F'(~,)
- F'(t~s)$1 = [fSL d3r'G°(r' r')(g'(qJ~(r'))-
g'(qJs(r'))~b(r'))[
< (sup fsL IGo(r, r')ld~r')(supLlg'( $~(r)) -- g'(t~s(r))l)l[t~[I.
(A16)
Let us now consider the ball S(0, R) C C(SL). From (A10) we have: Ig'(~b~) - g'(~s)[ = Kslf'(¢o + ~ ) - f'(~o + ¢2)[ × sup Ig'(t~(r)) - g'(¢s(r))[ r~SL
~< Ks sup_ If"(q~o + ~')1(11¢, - q'211). q~S(O, R)
Hence from (A16) we can write:
II(F'(q,,) - F'(q,=))q, II ~< Ks s u p \J~(f leo(r, r')l d~r'~follq,,/- t~2ll" [l~bll
(A17a)
86
J.T. DUNIEC AND S. W. THORNE fo =
sup
,/,E S(o, R)
If"(Oo+t~) I.
(A17b)
Thus the estimate of the constant L of the theorem can be expressed as:
L<~,~%sup(L, ICo(r,r')ld~r').
(A18)
Let K~= rain K~(Z) and K]= max K~)(z). Define the Green's functions G~(r,r') by 0
O
the equations: V2G~(r, r') - K~G~(r, r') = ~(r - r') i = 1, 2
(AI9)
with the boundary conditions G,(r, rs) = 0 i = 1, 2.
(A20)
By the strong maximum principle (Gilbarg and Trudinger, 1977) and the requirement that Gi(r,
r') r~r,)
exp ( ir - _K~lrr, I r'l) i = 1, 2,
-- --
(A21)
we find that G~(r, r') ~<0. The difference ~Go~ = G o - G~ satisfies the equation V28Gol
-
K~)SGoi =
( K oz-
K 2~ ) G , ,
(A22)
and hence by another application of the maximum principle we can show that G~ ~< Go ~< G3 ~<0. In view of these results, r~u~(fsL IG°(r' r')l d3r')~< rsup ( f ~ L - Gl(r, r')d3r').
(A23)
c l(r) = - fsL G l(r, r')d3r '
(A24)
The function:
is, by symmetry, independent of x and y and thus can be readily calculated by observing that it satisfies the equation: d2cl - K~cl = - 1, c,(O) = cl(d) = 0.
dz 2
(A25)
The actual form of the function c~(z) is: 1 shK~z _~ c,(z) = ~ ~ (exp (K,d) - 1) (exp (K,Z) -- 1).
(A26)
The maximum of this function is: 1 shKjd - (2(cosh K~d - 1)) ~j2 K~ shKld
max cl(z) = - ~
O~z~d
Using equations (A23, A24, A27) we can write the final estimate of the constant L:
(A27)
ELECTROSTATIC POTENTIALS IN MEMBRANE SYSTEMS
L<~sup([r6SL\JSL IG°(r' r')ld3r') ~<1~ ~ shK~d-(2(c°shK~d-1))'/2shKld
87 (A28)
For the application of the main theorem it remains to estimate the parameter rio, which in our case is: ri° = II-[Jz,=o "0~-'G°(r' ( r ' J d 2 r ~r'~) o z¢~s + f~.=a ~z,Go(r, r:)~b~ (r')d2r:[[.
(129)
By the reasoning similar to that applied above, we can show that f~5,O~(r, r')oz
<-"0-~TG°(r' r')oz
<~~ z 'G2(r' r') ~<0,
for z ' = 0
~z,G](r,r')~> 0~TGo(r,r')>~ 0-~G](r,r')~>0, for
(A30a)
z'=d.
(A30b)
With these results established we can perform the following estimates:
I - fz, oO-~,Go(r,r;)qJ,Ar'Jd2r'~+ f~,=do@,Go(r,r')O~Ar;)d2r'~ <~supl¢,,(rs)l[~.=o-O@,Gl(r,r')d2r',+fz.=a~z,G,(r,r's)d~r'].
(A31)
Define the functions:
u(z) = ~,
o-0@Gl(r, r')d2r:
(A32a)
v(z)= fz,=d~z,Gl(r,r')d2r's.
(A32b)
These functions obey the equations: d2u 2 (a) d ~ z - Klu = 0,
d 2v 2 (b) ~ - ~ - Klv = 0
(A33)
with the boundary conditions: (a)
u(d) = v(0) = 0,
(b) u(0) =
v(d) =
1.
(134)
The solutions to (A33, A34) are: (a) (b)
u(z) = [exp (K~(d - z)) - exp ( - K,(d - z))]/[exp (Kid) -- exp ( v(z) = [exp (K~z) - exp ( - K~z)]/[exp (K~d) - exp ( - K,d)].
rid)]
(135)
88
J.T. DUNIEC AND S. W. THORNE
Now combining the expressions (A35), (A31) and (A29), we obtain ~3° = SUPr SL I-- Jz'f =O ~0-'G°(r'azr's)O~s(r',)dZr's + f~,=a O~ G°(r' r's)+ls (r',)d2r'~
<~ sup [ @ , s ( r , ) [ s u p i u ( z ) + r s
v(z)l = IIo.ll,
(A36)
z
where we have introduced the 'surface' norm of qs~s as the supremum of its absolute values on the surface of the slab. Finally, grouping expressions (A36), (A28), (A18), (A15) and (A6), we obtain the estimate of the constant h from the theorem: ,,~ s h , , , d
h = FoLio ~< -~
K1
-
x/(2(cosh ,~,d - 1))4o,L,/,,2.,, shKld
(A37)
The above estimate concludes the derivation of the sufficient conditions for the convergence of the N e w t o n - K a n t o r o v i c h p r o c e d u r e for the Dirichlet problem of the nonlinear Poisson-Boltzmann equation in a slab. The convergence is assured if the constant h, bound by (A37), is less than 1]2.
APPENDIX
2
The function Q*-(z, a, "o) appearing in the solution of the linearized P o i s s o n - B o l t z m a n n equation is discussed in detail by Sauv6 and Ohki (1979). F o r completeness we quote the main points here. 1. The equation satisfied by Q*(z, a, 71) can be written as: 12+~z-d~2- Q+~ ~ ( 1
2 i f ( 1 + ~(l+rl/2)))dd~+Q+h(12)=0,2
(A38)
where h(12) = 6"00 + (4"q2 - 20-q - 8)11z - (16rl z + 4 n + 16)113 + (4n 2 - 20rl - 8)124 + 6-q125 [1 + (1 + 2rl)l) - (1 + 2rl)l) z - i23] 2 (A39) and 12
~/(exp (~0o( -+ a -+ d)) + •12) - ~/(n/2 + 1) ~/(exp (qSo( -!--a - d)) + hi2) ~ ~/(n/2 + 1) exp ( ~- K~/(1 + n/2)(z - (a + d))). (A40)
The above formulae should be read with either the upper sign or the lower sign. 2. The solution to (A38) can be obtained as a p o w e r series: Q-*(z, a, "q) = 2
aklIk,
(A41)
k=O
where ao = 1 and ak(k = 1, 2 . . . . ) are determined b y the following recursive equation:
ELECTROSTATIC POTENTIALS IN MEMBRANE SYSTEMS ak = S@k)[ak_6Sl(k - 6) + ak-sS2(k - 5) +
ak-4Sa(k
-
89
4)
+ ak-3S4(k - 3) + ak-2S3(k - 2) + ak-lS2(k - 1)].
(A42)
Si(k)(i = 1. . . . ,4) are in turn given by the following formulae: Si(k) = 0, for k < 0
(A43a)
S,(k) = k(k + 2Po)
(A43b)
S2(k) = k(k + 2Po)(47 + 2) + 6-0
(A43c)
S3(k) = k(k + 2Po)(4~/2- 1) - 8 - 20-9 + 4"o2
(A43d)
S4(k)
=
k(k + 2Po)( - 4 - 8"0 - 872) - 16 - 4 7 - 16-02,
(A43e)
where Po =
1+ ~(1 +
hi2).
(A44)
Equations (A38-A44) describe completely the functions Q+, for z > a + d, and Q - , for z>-a-d.
LITERATURE Barber, J., J. Mills and A. Love. 1977. "Electrical Diffuse Layers and Their Influence on Photosynthetic Processes". F E B S Lett. 74, 174-181. Barlow, C. A., Jr. 1970. "The Electrical Double Layer". In Physical Chemistry, I X A . Electrochemistry, Ed. H. Eyring, pp. 167-246. New York: Academic Press. Duniec, J. T., M. J. Sculley and S. W. Thorne. 1979. " A n Analysis of the Effect of Monoand Di-valent Cations on the Forces Between Charged Lipid Membranes with Special Reference to the Grana Thylakoids of Chloroplasts". J. theor. Biol. 79,473--484. and S. W. Thorne. 1979. " A n Explanation of the Proton Uptake of Chloroplast Membranes in Terms of Asymmetry of the Surface Charges". F E B S Lett. 105, 1--4. and - - . 1980. " A Theory of Charge Separation, Ion, Electron and Proton Transport in Photosynthetic Membranes Based on Asymmetry of Surface Charges". 2. theor. Biol. 85, 691-711. and - 1981. "Effects of Discrete Charges and Dielectric Properties of M e m b r a n e - W a t e r Interface on Electric Potentials inside M e m b r a n e s - - w i t h Reference to Charge Separation in Photosynthesis". F E B S Lett. 126, 1-4. Gilbarg, D. and N. S. Trudinger. 1977. Elliptic Partial Differential Equations of Second Order. Berlin: Springer Verlag. Krasnosel'ski, M. A., G. M. Vainikko, P. P. Zabreiko, Ya. B. Rutitski and V. Ya. Stetsenko. 1969. Approximate Solutions of Operator Equations. Moscow: Izd. Nauka. Levine, S., K. Robinson, G. M. Bell and J. Mingins. 1972. "The Discreteness-of-Charge Effect at Charged Aqueous Interfaces". J. electroanal. Chem. 38, 253-269. Loeb, A. L. 1951. " A n Interionic Attraction Theory Applied to the Diffuse Layer around Colloid Particles". J. Colloid Sci. 6, 75-91. Nakatani, H. Y., J. Barber and J. A. Forrester. 1978. "Surface Charges on Chloroplast Membranes as Studied by Particle Electrophoresis". Biochim. Biophys. A c t a 504, 215-225. Ninham, B. W. and V. A. Parsegian. 1971. "Electrostatic Potential Between Surfaces
-
-
-
-
90
J . T . D U N I E C A N D S. W. T H O R N E
Bearing Ionizable Groups in Ionic Equilibrium with Physiologic Saline Solution". J. theor. Biol. 31, 405-428. Sauvt, R. and S. Ohki. 1979. "Interactions of Divalent Cations with Negatively Charged Membrane Surfaces. I. Discrete Charge Potential". J. theor. Biol. 81, 157-179. Speck, S. H., S. Ferguson-Miller, N. Osheroff and E. Margoliash. 1979. "Definition of Cytochrome c Binding Domains by Chemical Modification: Kinetics of Reaction with Beef Mitochondrial Reductase and Functional Organization of the Respiratory Chain". Proc. natn. Acad. Sci. U.S.A. 76, 155-159. Waltz, D. 1979. "Thermodynamics of Oxidation-Reduction Reactions and Its Application to Bioenergetics". Biochim. Biophys. Acta 505, 279-353. RECEIVED 8-19-81 REVISED 1-14-82
0092-8240183]010069-22503.0010 Pergamon Press Ltd. © 1983 Society for Mathematical Biology
E L E C T R O S T A T I C P O T E N T I A L S IN MEMBRANE SYSTEMS • JACEKT. DUNIECand SIDNEYW. THORNE CSIRO Division of Plant Industry, P.O. Box 1600, Canberra 2601, Australia
A membrane with an arbitrary distribution of fixed charges inside and on its surfaces is considered. A procedure for calculating the local electrostatic potential at an arbitrary point of the system is described and its validity discussed. This procedure is based on the linearization of the 3-dimensional Poisson-Boltzmann equation around an exact 1dimensional solution.
1. Introduction. Several phenomena at a molecular level in biological membranes appear to be influenced by local electric fields. This is evident in particular in electron transfer reactions; for example, the redox reactions mediated by the cytochrome c depend significantly on the distribution of charges on this protein (Speck et al., 1979). Also, in photosynthetic membranes local transmembrane electric fields are likely to provide the necessary bias for the vectorial charge separation induced by light (Duniec and Thorne, 1979, 1980, 1981). Other examples are discussed in the review by Waltz (1979). It seems important, therefore, to be able to predict the magnitude of such fields in biomembranes as a function of charges inside and on the surface of the membrane as well as a function of ionic conditions in the surrounding solution. The simplest approach is based on the GouyChapman theory of electric potentials in electrolytes near charged surfaces of colloidal particles (Barlow, 1970). In this model one assumes that the membrane can be considered as an infinite flat layer sandwiched between the two half-spaces occupied by aqueous solutions. The surfaces of the membrane are charged with uniform surface densities which are representative of real densities of (discrete) fixed charges. This approach is adequate for the description of phenomena which depend on average electric potentials and fields. Examples can be: calculations of electrostatic repulsive force between the membranes (Duniec et al., 1979), local ionic concentrations near the membrane surface (Barber et al., 1977) and electrophoretic mobility of membranous particles (Nak69
70
J . T . DUNIEC AND S. W. THORNE
atani et al., 1978). However, for the processes at the molecular level the local electric fields and potentials are important, thus a model which accounts for the discrete nature of fixed charges and particular dielectric properties of the water-membrane interface is needed. Several approaches to build the model for the calculations of the discrete charge effects on the electric potentials in membranes have been made in the past [Sauv~ and Ohki (1979) and references therein]. The most recent method, developed by Sauv~ and Ohki (1979), allows for the calculation of electric potentials at the water-membrane interface for an arbitrary distribution of fixed charges on the membrane surface. Their procedure is based on the linearization of the Poisson-Boltzmann equation (describing electric potentials in electrolytes) around a certain 'smeared surface charge' solution. The Hankel transform technique is then applied to the resulting system of linear equations. A similar method was previously used in problems of adsorption of ions at the metal-water interface by Loeb (1951) and Levine et al. (1972). The work of Sauv6 and Ohki (1979) deals with a single water-membrane interface while the membrane is considered infinitely thick. Such an approximation is valid for calculations of local potentials due to charges near the surface of the membrane as there is only weak coupling across the membrane (Duniec and Thorne, 1981). The purpose of this paper is to provide a mathematical method for calculating local potentials anywhere in the membrane system. The fixed charges may be positioned on each side of the membrane as well as inside it. Such a method is needed for the problems related to the electron transfer reactions inside and across the membrane, where the internal charges and finite membrane thickness must be considered. An example of application of this method has already been published (Duniec and Thorne, 1981).
2. Mathematical Preliminaries. The fundamental physical law underlying the electrical phenomena in membrane systems is the Gauss' law of electrostatics: V. (e(r)V~(r)) = - 47rp(r),
(1)
where ((r), ~(r) and p(r) are the dielectric permittivity, electric potential and external charge density, respectively, at a point r = (x, y, z). The operator V = ((O/ax),(O/ay),(a/az)) and (.) denotes the scalar product. Throughout this work we shall be concerned with the system consisting of (see Figure 1): (1) th e aqueous region extending from z = - o o to z=-(a+d) andfrom z=a+dto z = + ~ ; ( 2 ) two polar regions: - a -
ELECTROSTATIC
POTENTIALS
ol° xy
I I I I
Wa~er
IN MEMBRANE
o o
SYSTEMS
71
I I
o
I o% I
o °
Membrane o
a 7-
Io
o 2d~
I o I I
----~ o
%
I
L I
o o
o
o
I I
o
I I
3
O
o
,:9
o
O
I
Figure 1. The coordinate system is chosen with the origin in the middle of the membrane and X, Y axes in its plane. The water region surrounds the membrane from both sides. Its interfaces with the membrane define the polar regions of thickness a. The dots in the membrane and the polar regions denote discrete charges.
d < z < -7 d and d < z < a + d; and (3) the hydrophobic portion of the membrane - d < z < d. The dielectric permittivities of each region will be assumed to be constant. Let us denote the dielectric constants as follows: ew for the water region, ep for the polar region and Em for the hydrophobic portion of the membrane. Let us first consider the aqueous region. We will restrict our attention to the case w h e n there are only salts of mono- and divalent cations dissolved in the water. All ionic species contribute to the volume charge density which, according to the Boltzmann's law, depend on the electric potential: p(r) = q[ - C - exp (qW/kT) + C + exp ( - q ~ / k T ) + 2C 2+ exp ( - 2q~/kT)],
(2) where q is the elementary charge, and C ÷, C - and C 2+ are the bulk concentrations of m o n o v a l e n t cations and anions and divalent cations respectively. F o l l o w i n g the notation of Ninham and Parsegian (1971), w e introduce the dimensionless electric potential tO= q'q/kT and the parameters 0 ~< rt ~< 1 and no such that: C 2+ = (1/2)~qno and C- = no. Since
72
J . T . D U N I E C A N D S. W. T H O R N E
the bulk solution must be neutral, i.e. 2 C 2+ + C + - C = 0, we can see that C + = (1 - "0)no. N o w f r o m (1) and (2) we obtain the following equation for tO(r): VztO = (/(2/2)[exp (tO) - (1 - "0) exp ( - tO) - "0 exp ( - 2tO)], or a + d < z < + o ~ .
-oo
(3)
The Debye's c o n s t a n t / ( is defined by /(2= 8.rrnoq2/ewkT. Equation (3) is a semilinear elliptic partial differential e q u a t i o n of the second order to which general t h e o r e m s of existence and u n i q u e n e s s of solutions to boundary value p r o b l e m s can be applied (Gilbarg and Trudinger, 1977). Practical m e t h o d s for solving such equations are, h o w e v e r , complicated. In A p p e n d i x 1 we shall discuss sufficient conditions for which the procedure based on the linearization of (3) a r o u n d an a priori known solution of the same equation is valid. Let tOo(Z) be a solution to (3), i n d e p e n d e n t of x and y, such that it vanishes at infinity. L e t
(4)
tO,,w(r) = tO(r) - tOo(Z).
Substituting (4) into (3), expanding (3) in terms of tO~,wand retaining only the leading terms, we obtain
2
2 '~.
v tO,, =/(o(Z).tO,.w,
(5)
where [in notation of Sauv6 and Ohki (1979)] Ko(Z) = (K2/2)[exp (tOo(Z) + (1 - n) exp ( - tOo(Z)) + 2"0 exp ( - 2tOo(Z))]
(z<-a-d
or z > a + d ) .
(6)
For notational c o n v e n i e n c e let us introduce the characteristic f u n c t i o n s for the polar regions and the interior of the m e m b r a n e :
X~(Z)
=fl'
t
X,,(z) X+ p(z)=
for - d - a < z < - d 0, elsewhere
=,[1, for d < z < d t 0, elsewhere
( 1,
ford
..
(7a)
(7b)
(7c)
E L E C T R O S T A T I C P O T E N T I A L S IN M E M B R A N E SYSTEMS
73
In addition, we shall use the superscript ( - ) for all the functions defined for z < d and the superscript ( + ) for those defined for z > d. As already implied in equations (4) and (7), the subscripts p, m and w will denote the quantities referring to the polar regions, membrane and water respectively. Let us assume that there are a total of N fixed charges at positions R~, and of valences v~ located in the membrane and the polar regions. Using equation (1) and the introduced notation, we can write the Poisson equation for the dimensionless potentials to~--in the polar regions and tom--in the membrane: N
V2to; = - ( 4 ¢ r q 2 / e p k T ) ~
v~8(r - r~)x~(z~)
(8)
v i S ( r - ri)xm(zj.
(9)
i=l
N
V2tom = - ( 4 7 r q 2 [ e m k T ) ~ i=l
The boundary conditions that are imposed on the potentials O?w, to~ and 0m can be derived from: (1) the continuity of the total electrostatic potential to(r) and the quantity e(dto/dz) (where E may be ew, Ep or Era) across the boundaries, and (2) the electroneutrality of the bulk solution which is equivalent to the requirement that the electrostatic potential vanished for z ~ +oo. These boundary conditions can be expressed as follows: lim O~,w(r)= 0
(lOa) +
too( +-- a +- d ) + to~,w(O, ++-a +- d ) = to;(O, +- a +- d )
Ew
+ Ew z=+-a±d
C
= Ep z=±
a±d
to;(p, _+ d) = to (p, -+ d)
(lOb)
z=±
a±d
(lOd) (lOe)
where p = (x, y). Linear equations (5, 8, 9) with the boundary conditions (10) can be solved by representing the potential in each region as the sum of contributions from individual charges (Sauv6 and Ohki, 1979) plus an additional 1-dimensional component due to a possible asymmetry between the potentials too(Z) and t0o(Z). Thus we seek the solution to equations (5, 8, 9) in the form:
74
J. T. DUNIEC AND S. W. THORNE N
O~,w(r) = ~bw(Z) + ~, ~w,~(r; r,)
(lla)
,=1
N
qJ;(r) = 6~(z) + Y. ~.,(r; r,)
(llb)
i=1
N
tOm(r) = q~,,(z) + ~ ~m,~(r; r,),
(llc)
i=l
where ~;,~, g~,~ and ~,~ are the individual charges contributions to the potentials 4'~.w, q'; and 4'm respectively, and qQ, 4~; and 4~,, are introduced to match the inhomogeneous boundary conditions (10b, c). The potentials ~,,, ~,~ and ~,~ obey the equations (i = 1, 2 , . . . , N): 2
m
2
±
Vr~,~(r; ri) = Ko(Z)~w,i(r; r,) 2
--+
--
4"a'q 2
(12a) +
V,~:p,,(r; r,) - - e p k T v~8(r - ri)x~(z~) Vr2G~.,(r; r,) -
4~rq 2 7~v,8(r-
r,)xm(z3.
(12b) (12c)
These potentials satisfy the following boundary conditions: ~,i(p, -+ a --- d; ri) = ~,~(p, -+ a -+ d; r~)
0G,,_ aG,,
ew0---7- - e p - ~ - , z = -+ a -+ d ~.,(p, ___d; ri) = ~m,,(p, -+d; ri)
(13a) (13b) (13c)
eo O~;'t aG,, oz = e m ~ , z = +- d
(13d)
lira ~w,~(r; ri) = 0.
(13e)
On the other hand, the potentials rbw(Z), 49;(z) and ~b,,(z) appearing on the r.h.s, of equation (11) obey the equations: ~ =
~:~(z)q~;
(14a)
ELECTROSTATIC POTENTIALS IN MEMBRANE SYSTEMS
d%; _d%m - ~ =0
75
(14b)
with the boundary conditions which depend on the choice of too(Z): &~( --- a + d) + too( ± a ± d) = ~b~(± a --- d) ew
+ew
=ep
, z=+a+d
&~(--- d) = 6•(± d) ep ~ z
(15a)
(15c)
~ =era dd4~ z' , z = ± d
(15d)
lira & ; ( z ) = 0.
(15e)
We have thus divided the linear system of equations (5, 8, 9) with the boundary conditions (10) into N independent systems (12, 13), each corresponding to an individual fixed charge in the polar region or inside the membrane, and one additional system (14, 15) whose solutions must be added to the sum of contributions of individual charges so that the total potential satisfies boundary conditions (10). Note that although boundary conditions (13) do not involve the function too(Z), the contributions of individual charges in each region depend on the choice of t0o(Z) through the function K~(z) [cf. eq (6)1 appearing in (12a). 3. C o n s t r u c t i o n o f t h e S o l u t i o n s to t h e L i n e a r i z e d P r o b l e m . We proceed now to construct the solutions to equations (12, 13) and (14, 15). Note that by symmetry, the potentials of individual charges depend on the absolute value of the difference p - p~, so that we can write:
•
~L(o, z;oi, z,) = ~;,,(R, z; z,)
(16a)
C A p , z; p,, z,) = C A R , z; z,)
(16b)
~,.,~(p, z; Oi, z,) = ~,.,~(R, z; zi),
(16c)
where
R=lp-p,I.
76
J.T.
DUNIEC
A N D S. W . T H O R N E
Accordingly, we can rewrite the operator V~ as: d 2
2
0 2
V~z = V~ + d-~, where VR = ~-~ + R -~ O0R '
(17)
and substitute equations (16, 17) into (12). In the case of (12a) this yields: /
I
J ~2x
n,~raz
I
(18)
= Ko(Z)~;.,(R, z ; z,).
"
Let us apply the Hankel transform defined as: dr
(19a)
f ( R ) = I o a f ( a ) S o ( a R ) da
(19b)
= Io Rf(R)Jo(aR)
to equation (18). After some algebra we arrive at: d. 2 ~+ 2 ~± d z 2 ¢~.i(a, z; zi) = (ct 2 + Ko(Z))~w,i(ol , z ; Zi).
(20)
Before applying a similar treatment to (12b, c) let us represent the functions ~.i and ~m.i as: 2 +
+
viq
+
~p.~(R, z; zi) = ~,,(R, z, z~) + epkTv,(R~ + (z - z~):) Xp(z~)
(21a)
2
+- z ')" ~m,i(R, z; z~) = ~,~(R, z, z~) + ~ , . k T v , ( R v~q ~ + (z - z~)~) X'n(
(21b)
Substituting (21) to (12b, c) and bearing in mind (16 and 17), we arrive at the following equations for ~.~ and ~m.~:
(
2
V2+
~p.~(R,z; z~) = 0
az /
"
Now the application of the Hankel transform to (22) yields:
(22a)
ELECTROSTATIC POTENTIALS IN MEMBRANE SYSTEMS
77
2~+
d ~19,i 2 "'+ dz 2 = a ~,~
(23a)
d2~,~,i 2dz = c~ ~,,,,i.
(23b)
Equations (20, 23) can be solved in a general form. In particular, the solution to (20) which vanishes at infinity can be written as [using the notation of Sauv6 and Ohki (1979)]:
~w,, = A~,iQ+-(z, a, "q) exp ( +- Paz),
(24)
where A~,, are certain constants to be found from the boundary conditions (13), and P is defined below: P = (1 + K2(1 + rl)/a2) '/2.
(25)
The functions Q-+(z, a, r/) (Sauv6 and Ohki, 1979) are described in Appendix 2. The general solutions to (23) can be written as:
~,, = A~,, exp ( - az) + B ~,, exp (az)
(26a)
~m,, = A,~,i exp ( - az) + B,,,., exp (az).
(26b)
In order to determine the 8 constants A~,,; A~.i; B~,i; Am,i; and Bin,i, we will derive a system of algebraic equations from the boundary conditions (13a-d). To this end, we use representations (16) and (21) in (13a-d) and then apply the Hankel transform (19) with respect to R = ]p -Oil. Using the well-known identity:
f R V/(R2 + 1(z -
zi) 2)
Jo(aR) d R
= 1 exp ( -
alz - z,I),
(27)
the above-described procedure yields: ~ , , ( a , ---a - d) = gp,i(a, -- a +- d) + /3, exp ( - ala + d -T-zil)X~(zi) ~ep (28a) d~w,, =
Ew d z
Ep
dd~/3,
exp ( - ala + d w z,l)Xp(z, )
(28b)
78
J . T . D U N I E C A N D S. W. T H O R N E
atz=
+_a+-d,
~;.~(a, +- d) + X~(Z,)a-~ exp ( - a]d -T z,]) = ~.,.,(a, +_ d) + X.,(z~) a ~ ~ exp ( - a ] d w- z~[) (28c)
e p ~ - - - / 3 i exp ( - aid -T- zil)X~(zi) = d~m.~ _
e,,--d~- + [3i exp ( - aid -~ zil)X,.(zi) (28d) at z = -+ d, where [3i = v~q2/kT. Note that equations (28) f o r m a s y s t e m of eight equations since all of the above formulae can be read with the u p p e r sign and with the lower sign. Let us introduce the eight-dimensional vector X consisting of the coefficients appearing in equations (24, 26): +
+
+
(X') r =(A~,,;A~,,;B~,,;Am,~;Bm,,;A~,,;Bp,~;Aw,~).
(29)
Let Y~ be the 8-dimensional vector w h o s e c o m p o n e n t s Y~, ] = l, 2 . . . . . 8 are defined as follows:
YI =
(13,1e,,a)x;(z,)
exp ( -
Y~2 = fl,X;(Z,) exp ( -
ala +
ala + d + z,I)
(30a)
d + z,I)
(30b)
Y~ = (fl,/a)(Xm(z,)/e,,, - X;(z,)/e~,) exp ( - aid + z,I)
(3Oc)
Y4 = (,Sffa)(X.,(z~) + X-p(zO) exp ( - aid + z,[)
(30d)
Y5 = ([3,/a)(xS(z,)/ep
- X . . ( z , ) / e . . ) e x p ( - aid - z,I)
(30e)
-b
Y6
=
(fl,la)(x,,(z,) + x,,,(z,))
exp
( -
a i d - z,I)
(30f)
Y~, = -(fl,/e~a)x;(z,) exp ( - ala + d - z,I)
(30g)
Y ; = ~,x+~(z,) exp ( - a l a + d - z,I).
(30h)
With the above defined notation, substitution of equations (24, 26) into
ELECTROSTATIC POTENTIALS IN MEMBRANE SYSTEMS
79
equation (28) results in a system of algebraic equations for the coefficients XI(j = 1 , 2 , . . , 8 ) which can be written in the following concise form: A X i = yi,
(31)
where the matrix A is defined by: p q 0 0 0 0 0 0
-
- -
-s
S
--1
r
A=
aeps t - ept 0 0 0 0
- a~ps t -1 ept -~ 0
-1
0
0
0
0
0
0
0
0
0
0
0
0
0
-
t
-
t -1
e,~t t -~
- e,~t -~ t
0
-- e,nt -1
emt
0 0
0 0
0 0
0 - t -~ ept ~ s -~ - epas -~
0 -t -- ept s %cts
0 0 0 p+ q+
(32) with the following notation: p-*(a) = 7 Q-*( -- a -- d, a, 7/) exp ( - P a ( a + d))
[ OQ*- 11
(33a)
-7- pQ*-( +_ a + d, a, ")In exp ( - P a ( a + d))
(33b) s ( a ) = exp ( a ( a + d))
(33c)
t ( a ) = exp ( a d ) .
(33d)
Upon solving (31) we obtain the coefficient vector Xi(a) [cf. (29)]. Thus we can find the Hankel transform of the potentials ~ , i ( a , z ; z~), ~,i(a, z; z~) and g,,,i(a, z; z~) from equations (24, 26). Subsequently, by inverting the Hankel transforms according to formula (19b) and using (21), we can find the ith charge contribution to the total electric potential in each region. In actual numerical calculations one needs to deal with the fact that for a = 0 the determinant of the matrix A [equation (32)] vanishes. However, the transforms ~;,~(a,z;zi), ~ . i ( a , z ; z O and ~,,,~(a, z; zi) are well defined for a = 0, which can be checked by examining the following limits: ~;,~(0, z; z~) = lim (A;,~(a)Q*-(z, a, "q) exp ( ~ P a z )
(34a)
80
J. T. DUNIEC AND S. W. THORNE
1
(;,i(~,z; zi) = lim Az,i(a)exp ( - az) + B z,i(a)exp (az) a-0
21.
+-J-exp
( - a i r - zil)x;(zi)]
E P ~
(34b)
{m,i(0, Z ;zi) = lim Am,i(a)exp ( - az) + Bm,i(a)exp (az) u -0
where A:,i, A;,i, B;,i, Am,i and BmViare the components of the vector x i ( a ) [equation (29)l. Performing the necessary algebra and using the notation
we obtain:
x
Q-(z7 O, ')
exp ( ~ ' ( 2+ a + d))
Q-( - a - d, 0, TJ)
(35c)
E L E C T R O S T A T I C P O T E N T I A L S IN M E M B R A N E S Y S T E M S
-aT~,-+~,,(0, ~l
81
z; z,) = (p -q ?P+_P +q-)(q-O + p-)
Q+(z, o, n)
x Q+(a + d, 0, a~) exp ( - K'(z - (a + d))).
(35g)
For the construction of the full solution to equations (8-10) we n e e d to find the potentials ~bw(Z), ~b~(z) and ~bm(z) appearing on the r.h.s, of equations (11) and satisfying equation (14) with the b o u n d a r y conditions (15). In general, these functions can be r e p r e s e n t e d as follows: ¢~,(z) = A~Q+-(z, 0, ~) exp ( ~ K'z) +
(36a)
&~(z) = A ~ z + B y
(36b)
(~m(Z) = A m z + Bin.
(36c)
Introducing the v e c t o r +
+
X ° = ( A : , A~, By, Am, Bin, Ap, Bp, Aw),
(37)
and using notation (33a, b), we can write b o u n d a r y conditions (15) as follows: --
O
0
0
p X 1 + X 2 ( a + d) - X 3 = - a -
0
0
--
(38a)
q X1 - epX2 = - b -
(38b)
- d X ° + X ° + dX°4 - X ° = 0
(38c)
epX ° - E~X°4 = 0
(38d)
dxo4+ X 5o - d X 6o -
(38e)
X ° = 0
e m X ° - epX ° = 0 +
0
p X8+(a+d)X°6+X +
0
0 6 =
(38f) °=a
+
+
(38g)
b ,
(38h)
a ± = tOo( -+ (a + d))
(39a)
q
X8
"~- E p X
w h e r e use was m a d e of the notation:
82
J . T . DUNIEC AND S. W. THORNE
b -+ =
ew
z = - (a + d).
(39b)
Once the solutions to (38) have been found, the full solution to the linearized problem can be obtained. As we have shown, this solution is determined by (a) the function too(Z), which is an exact solution to the 1-dimensional nonlinear Poisson-Boltzmann equation, and (b) the positions and valences of individual charges within the polar regions and the membrane. Our approach of approximating the solution to the 3-dimensional nonlinear equation (3) by a sum of an exact 1-dimensional solution to (3) and a solution to 3-dimensional linear equation (12a) is a first step in the Newton-Kantorovich method of constructing solutions to nonlinear problems (Krasnosel'ski e t a l . , 1969). In Appendix 1 we find the sufficient conditions under which this procedure as applied to the nonlinear Poisson-Boltzmann equation is convergent. Now we are concerned with finding an appropriate 'zeroth order' approximation t0o(z). An obvious choice is to let too(Z) be the solution to the 1-dimensional nonlinear Poisson-Boltzmann equation with the following boundary condition at z = - (a + d) (Sauv6 and Ohki, 1979): b ~- = ew
-
S
vi,
(40)
_+
where S is the area of the membrane, and I_+ are sets of indices denoting the charges in the r.h.s. ( + ) and the 1.h.s. ( - ) polar regions. Although this choice is not optimal, as there is some degree of coupling across the membrane, the possible discontinuity of the electric potential is rectified in the next step, i.e. the first-order approximation. We conclude this section by outlining a computational procedure which can be used in practical calculations. 1. Solve 1-dimensional Poisson-Boltzmann equation (3) with the boundary conditions (40). Obtain too(z). 2. Solve system of equations (38) using too(Z). Obtain ~bw(Z), $~(z) and qSm(z) according to (36, 37). 3. For each a needed to determine numerically the Hankel transform, and for each i, find the coefficient vector X i (29) from (31). Use an appropriate subroutine for solving a system of linear equations rather than explicit forms for X i, since these are very complex. 4. Using the coefficient vector X ~ (29) find the potentials ~,~, ~,i and ~,,.i by applying the inverse Hankel transform to equations (24, 26). For a ~ 1, use the asymptotic forms expressed by (35).
ELECTROSTATIC POTENTIALSIN MEMBRANESYSTEMS
83
5. The value of the total potential at a particular point r = (p, z) can be found by adding appropriate solutions obtained in steps 1-4. For example, the potential at the 1.h.s. polar region will be: N
4 (o, z) = 6 ; ( z ) + 5". ;(Io -o 1, z; z,) i=1
etc. 4. C o n c l u d i n g R e m a r k s . The purpose of this paper is to present a c o m p r e h e n s i v e description of the m e t h o d of calculating the local electric potentials in m e m b r a n e systems. We have avoided giving a specific example on purpose, since one application of the m e t h o d has already b e e n published (Duniec and Thorne, 1981). In this r e f e r e n c e the effects of salt conditions, positions of charges and dielectric properties of a w a t e r m e m b r a n e interface on the local potentials w e r e discussed. N e v e r t h l e s s it is worthwhile to m e n t i o n here the possible usefulness of the m e t h o d in our u n d e r s t a n d i n g of electron transfer processes in m e m b r a n e s . Electron transfer rates are affected by the electrostatic potential difference b e t w e e n the points w h e r e the donor and the acceptor reside. A change in this quantity m a y accelerate or slow d o w n the reactions. As we have shown r e c e n t l y (Duniec and Thorne, 1981), changes of the value of the dielectric permittivity of the polar regions as well as adsorption of ions m a y have a significant effect on the potential difference across the m e m b r a n e . It is possible, therefore, that the control of electron transfer rates in m e m b r a n e s m a y be exercised by such e n v i r o n m e n t a l changes. A n y relevant example would require a detailed k n o w l e d g e of the structure of the m e m b r a n e . As m o r e information on it b e c o m e s available, m a t h e m a t i c a l tools for quantifying electrical p h e n o m e n a at the molecular level will be needed. APPENDIX 1
The linearization of the Poisson-Boltzmann equation applied in this work is the first step in the Newton-Kantorovich method of solving nonlinear problems (Krasnosel'ski et al., 1969). Below we derive sufficient conditions for the applicability of this method to the nonlinear Poisson-Boltzmann equation and thus the validity of the linearization procedure. Consider the equation Fx = 0,
(A1)
where F is a Frechet ditterentiable nonlinear operator defined on a domain Dt of a Banach space B~ and transforming it into a Banach space B2. Let x* be a solution to (A1) and suppose that the nth approximation to x* has been found. Consider now the linear
84
J.T. DUNIEC AND S. W. THORNE
equation:
Fx. + F'(x.)(x - x.) = 0.
(A2)
If the linear operator [F'(x.)] ' exists then the solution to (A2) can be written as
x.+, = x. - [F'(x.)]-'Fx..
(A3)
The above formula defines the iterative Newton-Kantorovich method of solving (A1). The following theorem states sufficient conditions for the convergence of the above iterative process (Krasnosel'ski et al., 1969). THEOREM. Let F be determined in a certain ball S(xo, R) with the centre at Xo and the radius R, let its Frechet derivative F'(x) (treated as an operator defined on S(xo, R) and with values in the Banach space of continuous linear operators transforming B, into BE) satisfy the Lipschitz condition in S(xo, R): LIE'(X)
-
(A4)
f ' ( y ) l l ~ L I I x - YlI,
let Fo = I[f'(xo)-'ll, To = IIF'(Xo)-'FXoll
(AS)
ho = FoLno < 1/2.
(A6)
and
Then the sequence {x,},=o,, .... defined iteratively by (A3) converges to the solution x* E S(xo, R), where R > (1 - (1 - 2ho)'l:)rto/ho.
(A7)
We shall now apply this theorem to the Dirichlet problem of the nonlinear PoissonBoltzmann equation determined in a slab SL ={(x, y, z ) ; 0 < z
(A8)
with the boundary conditions: 0(rs) = Os(rs), where 0s is a bound, continuous function and r, = (x, y, 0) or (x, y, d). Let ~o(z) be a 1-dimensional solution to (A8), i.e. d20° = ~2f(Oo). dz 2
(A9)
Representing ~ as ~ = 6o + ~1, we can rewrite ( A l l ) as follows: V2~11-- K2(Z)I~I = K2 f(00 ~- ~1)-- f(o0) --
~Jl
--- g(O,),
(AI0)
where we have also defined the function g(t~,). If we let KZo(z)= K2f'(tko(Z)), then
g(o) = g'(O) = 0. Let Go(r, r') be the Green's function for the Dirichlet problem of the equation:
V20, - Ko2(Z)t~,= 0.
(All)
ELECTROSTATIC POTENTIALS IN MEMBRANE SYSTEMS
85
By means of Go(r, r') we can write (A10) in the integral form: tpl(r) = fSL d3r'G°(r' r')g(t~l(r')) - f~=o d2r5 0~Go(r, r's)t~ls(r's) + ~ =a dSrs 0 ~ G ( r , r's)t~,s(r's),
(A12)
where 6~s(rs) = ~ ~s(r~)- 60(0) for z = 0 b k , ( r ~ ) - ~o(d) for z = d" Equation (A12) defines the operator F acting on functions q'l: =- ~sL d3r'G°(r' r')g(t~(r')) + b(r) - t~x = 0.
(13)
The operator F acts in the Banach space of continuous functions C(SL) defined in the slab SL and equipped with the norm:
(A14)
I1~11= sup Iq,(r) I. O
The Frechet derivative of the operator F at some point $~E calculated from (A13):
C(SL)
can be readily
F'(t~l)t~ = fo
(A15a)
It is also apparent that F'(0)~b = - 6
(A15b)
r ° = II(F'(0)-')ll-- 1.
(MSc)
In order to estimate the Lipschitz constant of the operator F we perform the following:
I(F'(~,)
- F'(t~s)$1 = [fSL d3r'G°(r' r')(g'(qJ~(r'))-
g'(qJs(r'))~b(r'))[
< (sup fsL IGo(r, r')ld~r')(supLlg'( $~(r)) -- g'(t~s(r))l)l[t~[I.
(A16)
Let us now consider the ball S(0, R) C C(SL). From (A10) we have: Ig'(~b~) - g'(~s)[ = Kslf'(¢o + ~ ) - f'(~o + ¢2)[ × sup Ig'(t~(r)) - g'(¢s(r))[ r~SL
~< Ks sup_ If"(q~o + ~')1(11¢, - q'211). q~S(O, R)
Hence from (A16) we can write:
II(F'(q,,) - F'(q,=))q, II ~< Ks s u p \J~(f leo(r, r')l d~r'~follq,,/- t~2ll" [l~bll
(A17a)
86
J.T. DUNIEC AND S. W. THORNE fo =
sup
,/,E S(o, R)
If"(Oo+t~) I.
(A17b)
Thus the estimate of the constant L of the theorem can be expressed as:
L<~,~%sup(L, ICo(r,r')ld~r').
(A18)
Let K~= rain K~(Z) and K]= max K~)(z). Define the Green's functions G~(r,r') by 0
O
the equations: V2G~(r, r') - K~G~(r, r') = ~(r - r') i = 1, 2
(AI9)
with the boundary conditions G,(r, rs) = 0 i = 1, 2.
(A20)
By the strong maximum principle (Gilbarg and Trudinger, 1977) and the requirement that Gi(r,
r') r~r,)
exp ( ir - _K~lrr, I r'l) i = 1, 2,
-- --
(A21)
we find that G~(r, r') ~<0. The difference ~Go~ = G o - G~ satisfies the equation V28Gol
-
K~)SGoi =
( K oz-
K 2~ ) G , ,
(A22)
and hence by another application of the maximum principle we can show that G~ ~< Go ~< G3 ~<0. In view of these results, r~u~(fsL IG°(r' r')l d3r')~< rsup ( f ~ L - Gl(r, r')d3r').
(A23)
c l(r) = - fsL G l(r, r')d3r '
(A24)
The function:
is, by symmetry, independent of x and y and thus can be readily calculated by observing that it satisfies the equation: d2cl - K~cl = - 1, c,(O) = cl(d) = 0.
dz 2
(A25)
The actual form of the function c~(z) is: 1 shK~z _~ c,(z) = ~ ~ (exp (K,d) - 1) (exp (K,Z) -- 1).
(A26)
The maximum of this function is: 1 shKjd - (2(cosh K~d - 1)) ~j2 K~ shKld
max cl(z) = - ~
O~z~d
Using equations (A23, A24, A27) we can write the final estimate of the constant L:
(A27)
ELECTROSTATIC POTENTIALS IN MEMBRANE SYSTEMS
L<~sup([r6SL\JSL IG°(r' r')ld3r') ~<1~ ~ shK~d-(2(c°shK~d-1))'/2shKld
87 (A28)
For the application of the main theorem it remains to estimate the parameter rio, which in our case is: ri° = II-[Jz,=o "0~-'G°(r' ( r ' J d 2 r ~r'~) o z¢~s + f~.=a ~z,Go(r, r:)~b~ (r')d2r:[[.
(129)
By the reasoning similar to that applied above, we can show that f~5,O~(r, r')oz
<-"0-~TG°(r' r')oz
<~~ z 'G2(r' r') ~<0,
for z ' = 0
~z,G](r,r')~> 0~TGo(r,r')>~ 0-~G](r,r')~>0, for
(A30a)
z'=d.
(A30b)
With these results established we can perform the following estimates:
I - fz, oO-~,Go(r,r;)qJ,Ar'Jd2r'~+ f~,=do@,Go(r,r')O~Ar;)d2r'~ <~supl¢,,(rs)l[~.=o-O@,Gl(r,r')d2r',+fz.=a~z,G,(r,r's)d~r'].
(A31)
Define the functions:
u(z) = ~,
o-0@Gl(r, r')d2r:
(A32a)
v(z)= fz,=d~z,Gl(r,r')d2r's.
(A32b)
These functions obey the equations: d2u 2 (a) d ~ z - Klu = 0,
d 2v 2 (b) ~ - ~ - Klv = 0
(A33)
with the boundary conditions: (a)
u(d) = v(0) = 0,
(b) u(0) =
v(d) =
1.
(134)
The solutions to (A33, A34) are: (a) (b)
u(z) = [exp (K~(d - z)) - exp ( - K,(d - z))]/[exp (Kid) -- exp ( v(z) = [exp (K~z) - exp ( - K~z)]/[exp (K~d) - exp ( - K,d)].
rid)]
(135)
88
J.T. DUNIEC AND S. W. THORNE
Now combining the expressions (A35), (A31) and (A29), we obtain ~3° = SUPr SL I-- Jz'f =O ~0-'G°(r'azr's)O~s(r',)dZr's + f~,=a O~ G°(r' r's)+ls (r',)d2r'~
<~ sup [ @ , s ( r , ) [ s u p i u ( z ) + r s
v(z)l = IIo.ll,
(A36)
z
where we have introduced the 'surface' norm of qs~s as the supremum of its absolute values on the surface of the slab. Finally, grouping expressions (A36), (A28), (A18), (A15) and (A6), we obtain the estimate of the constant h from the theorem: ,,~ s h , , , d
h = FoLio ~< -~
K1
-
x/(2(cosh ,~,d - 1))4o,L,/,,2.,, shKld
(A37)
The above estimate concludes the derivation of the sufficient conditions for the convergence of the N e w t o n - K a n t o r o v i c h p r o c e d u r e for the Dirichlet problem of the nonlinear Poisson-Boltzmann equation in a slab. The convergence is assured if the constant h, bound by (A37), is less than 1]2.
APPENDIX
2
The function Q*-(z, a, "o) appearing in the solution of the linearized P o i s s o n - B o l t z m a n n equation is discussed in detail by Sauv6 and Ohki (1979). F o r completeness we quote the main points here. 1. The equation satisfied by Q*(z, a, 71) can be written as: 12+~z-d~2- Q+~ ~ ( 1
2 i f ( 1 + ~(l+rl/2)))dd~+Q+h(12)=0,2
(A38)
where h(12) = 6"00 + (4"q2 - 20-q - 8)11z - (16rl z + 4 n + 16)113 + (4n 2 - 20rl - 8)124 + 6-q125 [1 + (1 + 2rl)l) - (1 + 2rl)l) z - i23] 2 (A39) and 12
~/(exp (~0o( -+ a -+ d)) + •12) - ~/(n/2 + 1) ~/(exp (qSo( -!--a - d)) + hi2) ~ ~/(n/2 + 1) exp ( ~- K~/(1 + n/2)(z - (a + d))). (A40)
The above formulae should be read with either the upper sign or the lower sign. 2. The solution to (A38) can be obtained as a p o w e r series: Q-*(z, a, "q) = 2
aklIk,
(A41)
k=O
where ao = 1 and ak(k = 1, 2 . . . . ) are determined b y the following recursive equation:
ELECTROSTATIC POTENTIALS IN MEMBRANE SYSTEMS ak = S@k)[ak_6Sl(k - 6) + ak-sS2(k - 5) +
ak-4Sa(k
-
89
4)
+ ak-3S4(k - 3) + ak-2S3(k - 2) + ak-lS2(k - 1)].
(A42)
Si(k)(i = 1. . . . ,4) are in turn given by the following formulae: Si(k) = 0, for k < 0
(A43a)
S,(k) = k(k + 2Po)
(A43b)
S2(k) = k(k + 2Po)(47 + 2) + 6-0
(A43c)
S3(k) = k(k + 2Po)(4~/2- 1) - 8 - 20-9 + 4"o2
(A43d)
S4(k)
=
k(k + 2Po)( - 4 - 8"0 - 872) - 16 - 4 7 - 16-02,
(A43e)
where Po =
1+ ~(1 +
hi2).
(A44)
Equations (A38-A44) describe completely the functions Q+, for z > a + d, and Q - , for z>-a-d.
LITERATURE Barber, J., J. Mills and A. Love. 1977. "Electrical Diffuse Layers and Their Influence on Photosynthetic Processes". F E B S Lett. 74, 174-181. Barlow, C. A., Jr. 1970. "The Electrical Double Layer". In Physical Chemistry, I X A . Electrochemistry, Ed. H. Eyring, pp. 167-246. New York: Academic Press. Duniec, J. T., M. J. Sculley and S. W. Thorne. 1979. " A n Analysis of the Effect of Monoand Di-valent Cations on the Forces Between Charged Lipid Membranes with Special Reference to the Grana Thylakoids of Chloroplasts". J. theor. Biol. 79,473--484. and S. W. Thorne. 1979. " A n Explanation of the Proton Uptake of Chloroplast Membranes in Terms of Asymmetry of the Surface Charges". F E B S Lett. 105, 1--4. and - - . 1980. " A Theory of Charge Separation, Ion, Electron and Proton Transport in Photosynthetic Membranes Based on Asymmetry of Surface Charges". 2. theor. Biol. 85, 691-711. and - 1981. "Effects of Discrete Charges and Dielectric Properties of M e m b r a n e - W a t e r Interface on Electric Potentials inside M e m b r a n e s - - w i t h Reference to Charge Separation in Photosynthesis". F E B S Lett. 126, 1-4. Gilbarg, D. and N. S. Trudinger. 1977. Elliptic Partial Differential Equations of Second Order. Berlin: Springer Verlag. Krasnosel'ski, M. A., G. M. Vainikko, P. P. Zabreiko, Ya. B. Rutitski and V. Ya. Stetsenko. 1969. Approximate Solutions of Operator Equations. Moscow: Izd. Nauka. Levine, S., K. Robinson, G. M. Bell and J. Mingins. 1972. "The Discreteness-of-Charge Effect at Charged Aqueous Interfaces". J. electroanal. Chem. 38, 253-269. Loeb, A. L. 1951. " A n Interionic Attraction Theory Applied to the Diffuse Layer around Colloid Particles". J. Colloid Sci. 6, 75-91. Nakatani, H. Y., J. Barber and J. A. Forrester. 1978. "Surface Charges on Chloroplast Membranes as Studied by Particle Electrophoresis". Biochim. Biophys. A c t a 504, 215-225. Ninham, B. W. and V. A. Parsegian. 1971. "Electrostatic Potential Between Surfaces
-
-
-
-
90
J . T . D U N I E C A N D S. W. T H O R N E
Bearing Ionizable Groups in Ionic Equilibrium with Physiologic Saline Solution". J. theor. Biol. 31, 405-428. Sauvt, R. and S. Ohki. 1979. "Interactions of Divalent Cations with Negatively Charged Membrane Surfaces. I. Discrete Charge Potential". J. theor. Biol. 81, 157-179. Speck, S. H., S. Ferguson-Miller, N. Osheroff and E. Margoliash. 1979. "Definition of Cytochrome c Binding Domains by Chemical Modification: Kinetics of Reaction with Beef Mitochondrial Reductase and Functional Organization of the Respiratory Chain". Proc. natn. Acad. Sci. U.S.A. 76, 155-159. Waltz, D. 1979. "Thermodynamics of Oxidation-Reduction Reactions and Its Application to Bioenergetics". Biochim. Biophys. Acta 505, 279-353. RECEIVED 8-19-81 REVISED 1-14-82