Abstract:

Milk is a complex liquid, which contains many different species, for example proteins, fat, minerals etc. It is the primary source of
nutrition for young mammals before they are able to digest other types of food. The proteins can be divided into two groups: caseins and
whey. Whey proteins are about 20 wt % of the total protein amount in milk, whereas the caseins corresponds 80 wt % of the total protein
content in milk. The largest structures in the fluid portion of the milk are “casein micelles” which are aggregates of several thousands of
protein molecules. The micelle is considered to be spherical and the diameter is in the micrometer size.

The caseins can be divided in four types:
α
s1- ‐casein,
α
s2- ‐casein,
α
- ‐casein, and
α
- ‐ casein. ß-
‐casein is one of the most abundant
caseins and it also self-
‐assembles to larger aggregates. In this thesis we have used a simple model to try to capture how electrostatic
interactions affects the structure of
α
- ‐casein micelles in milk. The micelles have been modeled as hard spheres, with a central
net charge of -
‐ 140e, and a radius of 75 Å. These parameters have been taken from experimental data published in the literature. The
structure of the solution has been studied by comparing the radial distribution functions for different solution conditions, such as the salt
concentration and valency, pH, and the temperature.
Popular scientific description

Milk is a complex liquid, which contains many different species, for example proteins, fat, minerals etc. It is the primary source of
nutrition for young mammals before they are able to digest other types of food. The proteins can be divided into two groups: caseins and
whey. Whey proteins are about 20 wt % of the total protein amount in milk, whereas the caseins corresponds 80 wt % of the total protein
content in milk. The largest structures in the fluid portion of the milk are “casein micelles” which are aggregates of several thousands of
protein molecules. The micelle is considered to be spherical and the diameter is in the micrometer size.

In this thesis we have used a simple model and computer simulations to try to capture how electrostatic interactions affects the
structure of
α
- ‐casein micelles in milk. The micelles have been modeled as hard spheres, with a central net charge of -
‐140e, and a
radius of 75 Å. The volume fraction of was set to 5%, which is the actual volume fraction in the real product. The structure of the
solution has been studied by comparing the radial distribution functions for different solution conditions, such as the salt concentration
and valency, pH, and the temperature.

It was noticed that due to the fact that the micellar charge is very large, the electrostatic repulsive interaction dominates, and the mean
distance between the micelles are almost always obtained. Moreover, an increase of the temperature does not affect the structure at
all i.e. the entropic contribution due to increased temperature can be neglected in comparison with electrostatic repulsion between the
micelles. Also, there might also be an influence of the electric permittivity since it was kept constant during the simulations,When the salt concentration was increased to 80 mM, which corresponds to the ionic strength in milk, the structure of the
α
- ‐casein
micelles resembles the structure of an ideal gas i.e. the electrostatic repulsive interactions are screened.

Introduction:

Milk is as ancient as human beings themselves. Historians believe
that humans started to drink milk over 10,000 years ago, along
with the start of the domestication of animals. Milk is the normal
product of the mammary glands of female mammals. Its purpose
is primarily to meet the nutritional requirements of the neonate.
It contains more than 100 substances, and it is the most complete
and complex known natural food. Milk is perishable because it is
a very good environment for the growth of the microorganisms
that caused the FDA to warn people from drinking raw milk in
March 2007. They wrote: Pasteurized milk does safe live! [41]
The pasteurization of milk is heating the milk in a certain way to
save the nutritional value of the milk while destroying the harmful
bacteria.

The tropical countries are not considered consumers of milk,
because of the high temperature and the lack of refrigeration, which
makes the storage of milk difficult. The north part of the world,
especially northern Europe and North America, is considered a big
consumer of milk and milk products. Milk has traditionally been
preserved by converting it to more stable products like making
cheese and yogurt by milk fermentation. [33]

1:1 Milk processing :

1:1:1 Pasteurization :

Because of the perishability of milk, milk storage has been an
important process. The most widely used process, for preserving
fresh milk, is pasteurization, which was named after the French
microbiologist Louis Pasteur in 1864.

The most popular methods of pasteurization are:

1) HTST: high temperature short time pasteurization, which
requires heating milk to 720 C for 15 seconds, and can extend the
shelf life to a week.
2) UHT: ultra high temperature pasteurization, which requires
heating the liquid to 1380 C for 2 seconds, and can extend the
shelf life for weeks. It causes some changes of the structure of the
milk protein.[2]

1:1:2 Fermentation :

Other processes, that are popular for preserving milk, involve
preserving milk components like fat and protein. These processes
are modifying the environment of the milk and extend the lifetime
of the products. Adding lactic acid bacteria to the milk converts
lactose (milk sugar) to lactic acid, which leads to the
dissociation of the casein micelles. After that, the viscous gel
solution has the sour taste of yogurt. After adding the lactic acid
bacteria to the milk, the lactose is converted to lactic acid, which
decreases the pH of the milk. Rennet is added to convert the casein
micelle into a stable protein that does not dissociate in water.

The content of the milk differs according the mammals it comes
from as shown in table 1. The components also depend on the
nutrition and the season. Cow’s milk for example contains about
87% water. [34]

1:2 Milk components :

1.2.1 Milk fat :

Fat represents about 3.5-
‐6.5 % of milk, and it is the most variable
component. It is secreted as a fat globule surrounded by a lipid
bilayer membrane. When raw milk is centrifuged, the fat moves
to the top and forms a cream layer, because the lipid is lower in
buoyant density, see Figure 1

Figure 1: shows the formation of two layers upon centrifugation of
the milk. The upper layer is the lipid and it forms a creamy layer,
and the layer at the bottom is formed from skim milk (milk without
fat). This is because the lipid layer has a low buoyant density.

The cream contains some proteins (Mucin, Xanthine oxidase and
butyrophilin), which are carried by the fat globules. It has been
suggested that the fat globules are accumulated between the two
halves of the milk lipid globule membrane (MLGM), which is part
bilayer and part protein (see figure 2). Milk fat globules (MFG),
whose diameters range from 0.2 μm to 20 μm, with a mean
diameter of 4μm (Mulder and Walstra, 1974), are essentially from
triglycerides, which are esters derived from glycerol combined
with three fatty acids. It has been found that the triglyceride forms
98% of the fat content, and the phospholipid forms 1% of the fat
content.

Over 400 fatty acids have been identified in milk fat. The most
common type of fatty acids is triglycerol or triglyceride, see
figure 3, which is composed of three fatty acids.

Milk fat melts over a wide range of temperatures from -
‐400C to
400C, and it can be degraded by excess heat or light, and actions
of some enzymes cause it to be inactivated.

1.2.2 Lactose :

Lactose is the major carbohydrate in milk, and it is a disaccha
-
ride derived from galactose and glucose as seen in Figure 4. It is
found dissolved in the milk serum in two forms: alpha-
‐anomer
and beta-
‐anomer. The ratio between the two forms changes with
temperature. The beta-
‐anomer is sweeter and more soluble
than the alpha-
‐anomer. Short pasteurization time does not have a
significant effect on the glucose, but (UHT pasteurization) moti
-
vates the Brownian reaction between the Lactose and the protein,
leading to undesirable flavor and color

1.2.3 Minerals :

The majority of minerals found in milk are Ca, P, Mg, K and
Zn. Calcium and phosphorus combine together to form Calcium
Phosphate salt and it is called Colloidal Calcium Phosphate
(CCP). It was suggested to play an important rule in binding the
casein micelle by immigrating in or out the micelle with changing
temperature.

1.2.4 Vitamins

Milk is not considering a major source of vitamins in a diet. There
are two groups of vitamins that exist in milk according to their
solubility:
•
The water-
‐soluble vitamins group, including the
vitamins (B1, B2, B3, B5, B6, B12, C).
•
The fat-
‐soluble vitamins group, including the vitamins
(A, D, E, K).
HTST does not affect the vitamins, but UHT decreases
the water-
‐soluble vitamins. Exposure to light affects the vitamin
A content.

1.2.5 Milk proteins :

The proteins in milk are unique, and make milk an important
source of human nutrition. Milk proteins are digestible by almost
everyone in the intestine. Generally, there are many types of
proteins according to the sequences of the amino acids. There are
20 amino acids, which form the different protein types where 9
of them are essential for diet and all of which are present in milk.
The amino acids are bonded to each other by polypeptide bonds,
see Figure 5. The N-
‐ terminal (terminus) is the amine part and
the C-
‐terminal (terminus) is the carboxylic acid part. They
carry negative charge (the N terminal) and positive charge (C-
‐
terminal), and they neutralize each other. If the hydrocarbon
side chain can be titrated, the charge of the protein depends on the
pH, which will be discussed in detail later.

A denatured protein is the protein unfolded from its native state.
The denaturation can be useful in the industry or digestion as it
gives the enzymes access to the protein chains.

Bovine milk proteins are the most common proteins. Caseins
resemble about 80% of all the protein and it is present in milk
of all species. All other proteins are grouped together in a group
called whey proteins. The major whey protein in cow milk is
alpha-
‐lactalbumin and beta-
‐lactoglobulin. The only known
milk allergy is caused by beta-
‐lactoglobulin, because of its
poor digestibility. The whey protein has also an important use in
the industry as it is used in binding meat and water in meat and
sausage products, and it will not form a gel upon denaturation and
acidification.

There are four types of caseins in milk, each having appropriate
amino acids that make the milk essential for growth and
development of the young. They they are also phosphorylated to
various degrees, which plays an important role in the stabilization
of the casein itself.

A casein micelle is composed of several similar proteins,
forming a multi-
‐ molecular granular structure. The casein micelle
contains water plus salts (mainly from Calcium and Phosphorus).
The casein micelle plays an important role in many milk product
industries such as in cheese making. The casein can be separated
from whey by centrifugation where the casein is pelleting and the
whey is in supernatant, or it can be precipitated by adding acid.

The average diameter of a casein micelle is around 120 nm to
180 nm, and it is about 1/50 of the fat globule size. The micelle
is composed of proteins with salt and the salt is primarily nano-
‐
clusters of Calcium and Phosphorus. Researchers think that the
Calcium and Phosphorus play a role in connecting the proteins
because by removing the salts, the micelles dissociate into small
parts called sub-
‐ micelles. The micelles were found to have
essentially 4 proteins;
α
- ‐s1,
α
- ‐s2, ß and χ-
‐casein in the ratio
4:1:4:1. [2] The isoelectric point of the casein is at pH

4.6 i.e. the pH at which the net charge is zero. At pH 6.6, the
nano-
‐clusters, which are mainly Calcium and Phosphorus salts,
are dispersed as small particles with a radius of 1.72 nm-
‐ 2.27
nm. There are many articles about the internal structure, and all of
them agree that there is a micelle formed by the proteins and the
colloidal calcium phosphate (CCP). [27]–[30]

Generally, the studies that use light and x-
‐ray and neutron
scattering show that the micelle is formed from sub-
‐units (called
sub- ‐micelles) with the help of CCP. CCP acts like cement, as
it connects the sub-
‐micelles. The studies that use electron
microscopes show that there are no sub-
‐micelles, but there is a
protein matrix and channels in the inside of the micelle. However,
some recent studies, using high-
‐resolution microscopes and
carefully made samples show that the sub-
‐micelle model
probably exists. Horne (2006) has given an explanation, to the
sub- ‐micelle model, as is shown in Figure 6, and illustrates that
the whole micelle consists of small spherical sub-
‐units (sub-
‐
micelles) connected with the colloidal calcium phosphate. The
sub- ‐micelles on the outer surface of the micelle are rich with k-
‐
casein, having protruding hairs of negative charges which induce
a steric repulsion between micelles, thus preventing them to
coagulate and hence stabilizes the milk.

The Holt model explains how the nano-
‐cluster of the CCP acts
like cement, and connect the proteins to form the big micelle as is
shown in Figure 7. Both alpha and beta casein are attached to the cluster. The difference between
the proteins enables them to make a bridge with the cluster.

Old electron microscope studies say that there is no inhomogeneity
in the micelle more than a few nanometers, which means that
all the components form one part, and that there are no small
molecules looking like the big micelle. However, the previous
electron microscope photos contain artifacts, since the samples
were supposed to be fixed by drying or freezing, with a possible
loss of protein flexibility.

Old electron microscope studies say that there is no inhomogeneity
in the micelle more than a few nanometers, which means that
all the components form one part, and that there are no small
molecules looking like the big micelle. However, the previous
electron microscope photos contain artifacts, since the samples
were supposed to be fixed by drying or freezing, with a possible
loss of protein flexibility.

In general all studies confirm that there is a casein micelle,
composed of four types of protein, namely alpha s1, alpha
s2, beta and kappa casein. The protein micelle is similar to the
surfactant micelle with both a hydrophobic part and a hydrophilic
one, motivating the micelle formation by their behavior towards
the water phase. The polar part, which is presented by the N
and C terminals, tends to face the water. The hydrophobic part
is presented by two portions, the first of which consisting of non
polar side chains and the second is the repeated sequence of P, F

(Phenylalanine and Proline). It has been found that P, F sequences
have less polarity than other sequences.Milk can be considered as a colloidal solution that contains:
•
Solution from milk sugar in water.
•
Fat emulsions (oil in water).
•
Proteins suspended in water.

Because of the great role played by the casein micelle, a simple
model has been used to simulate the role of electrostatic interactions
between ß-
‐casein micellar structures in milk with emphasize on
the following:1)
The effect of the electrostatic environment on the
micelle-
‐micelle interaction.
2)
The effect of temperature of the micelle-
‐micelle
interaction.
3)
The effect of changing pH on the protein charge on
the micelle-
‐micelle interaction.
2.
Theoretical Background,
2.1
Electrostatic interactions,

Milk contains approximately 87% percent water. It is extremely
time consuming to consider the effects of all the individual water
molecules in computer simulations, and therefore we will only
consider the global effect of water, by letting the electrostatic
interaction energy (uel) between charges, be screened by the
dielectric permittivity of water, denoted Ɛr.

2.1, 2.2

The electrostatic energy (푢!”)
depends
on the electrostatic
potential
(Φ!”). The electrostatic potential depends on the charges found
in the system. The electrostatic potential in vacuum at a point,
separated by a distance rj from a point charge j, is given by:

2.3

where Ɛ! is the dielectric permittivity of vacuum. The
electrostatic interaction energy between two particles i and j (in
vacuum) can be obtained by:

2.4

This equation represents Coulomb’s law. Coulomb’s equation for
electrostatic interaction in water is given by:

2.5

As mentioned above, the introduction of the dielectric permittivity
coefficient of water Ɛ!, accounts for the screening of the
electrostatic interaction between charges when treating water
molecules implicitly as a structureless continuum.

The effect of the concentration, the valence of the ions, and
the temperature on the system, can be studied using the model
(explained above). The Debye length can be considered as a
resisting factor of the electrostatic force and is defined as
length
. The Debye length can be depicted as the radius of the
circle with a
point charge oriented in its origin, and the circumference of this
circle is the farthest position that can be influenced by this point
charge

2.6

The Debye length affects the electrostatic potential, and the
electrostatic potential from particle j, which affect particle i,
becomes,

2.7

which is the Debye Hückel expression for the electrostatic potential.
Here Rj is the radius of the particle, and rij is the distance between
i and j. It is obvious that the Debye screening length decrease by
the increasing of the concentration of salt and the charge, but also
while increasing temperature.

In our model, the radial distribution function between the micelles
has been calculated numerically to show the electrostatic behavior
of the model when changing the environment.

2.2 Protein charge :

The N and C terminus, and several side chains of a peptide group
can be charged. The protein charge can be estimated by calculating
the sum of the amino acids charges. The terminals of the amino
acids, which are an amine and a carboxylic acid group, compensate
each other, so they can be considered as neutral. The side chain
can be charged and it can be dissociated at a suitable pH. The
acidic groups (AH) dissociate into a proton (H+) and the conjugate
negative ion (A-
‐):

The basic groups (BoH) dissociate into a hydroxyl group (OH-
‐)
and a positive ion
(B+):

The equations (1.1) and (1.2) can be combined as

This equation can be generalized for any titration process as B
is the titrand with the charge ne against a titrant with a charge
me, resulting in ABP where p= n+ m, where (e) is the elementary
charge.

For this reaction the dissociation constant is:

which is equivalent to:

By assuming that the activity coefficient is unity, if An is a proton
and Kd is referred to acid dissociation constant Ka, pA = pH

x
!
!
!
is obtained by the same steps
.

When pH is equal to pKd, the exponents in the mole fractions
equations will be zero, which leads to be ½ for both the titrant and
the titrand. When pH pKd, there are an increase of the number
of titrant sites. For pH >pKd the number of titrant sites decreases.

It is possible afterward to express the total charge by mole
fractions at constant pH by summing x!! over all the basic groups,
substracting them from the sum of
x!”! over all the acidic groups.

3. Model and Method :

3:1 Model :

A simple model system has been designed to study the micelle-
‐
micelle interactions in milk. The charge of the proteins which
build up the micelle, can be estimated from the following steps:
1.
Finding the sequence of amino acids of the protein. [36]
2.
Picking up the polar chains and knowing if they are
positive or negative according the pH. [37]
3.
Substituting in equation (2.22)
At a temperature of 250C and pH=6.7 (conditions of fresh milk)
the net charge of a micelle protein is about -
‐7 elementary charges.
In our model, we treated the micelles as hard spheres with a radius
of 75 Å, and we neglected any titration effects during the coarse of
the simulations. Practical experiments have been done to calculate
the micellization number of the ß casein.[38] It is found to be
around 20, i.e., there are 20 proteins/micelle, and the total charge
would thus be around (-
‐140) e/micelle.

Figure 9:
Schematic draw of the reference model
consisting of 21 micelles and 2940 counterions. All
species are modeled as hard spheres, and the water is
treated by applying a dielectric continuum.

In order to make the system overall electroneutral, 140
counterions/micelle were added to the system. The ions were
modeled in the same way as the micelles, i.e. as point charges,
enclosed in impenetrable spheres, but with a diameter of 4 Å. The
smallest investigated system consisted of 21 micelles and 2940
counterions, enclosed in a cubic simulation box with side lengths
equal to 906 Å. The dimensions of the box were chosen as to give
a micelle volume fraction of about 5 %, a value corresponding to
the conditions in milk. The volume fraction was simply calculated
as the total volume of the spherical model micelles, divided by
the box volume. A snapshot of the model system is provided in
Figure 8.

3:2 Method :

Fortunately, most of the computer simulations are based on the
assumption that classical mechanics can be used to describe
the motion of atoms and molecules. The fact that the quantum
system can be found in different states is needed. Some important
definitions with respect to computer simulations are:

International Journal of Biotechnology and Bioengineering
Volume 3 Issue 6, July 2017
209
assumption that classical mechanics can be used to describe
the motion of atoms and molecules. The fact that the quantum
system can be found in different states is needed. Some important
definitions with respect to computer simulations are:
Ergodic system: when the time average of the sequence of events
is equal to the ensemble average.
The difference between analytical and numerical calculations: the
analytical calculations are exact and can be done by using pen and
paper and the solutions are exact. The numerical calculations are
more complicated and can be done by using computer simulations,
accompanied with an inevitable statistical noise.

Markov chain: a system is said to undergo Markov chains when
it undergoes a chain of transitions among finite possible states
and the next state only depends on the current state and not the
sequence of the state.
Importance sampling: A technique used for estimating a particular
distribution for a sample generated by random distribution. For
example, two distributions as shown in Figure 10 where f(x) is the
random distribution, and h(x) is the particular one ’’certain one’’.

Importance weights: It is a function used to measure the error
between the two distributions,
=
3.1
hence the rejection rate depends on the importance weights. We
need a new distribution related to h(x) only. By using the importance
weight technique, approximation of the new distribution can be
obtained.

3.1

3.2

where E is called the Monte Carlo estimator. The variance can be
calculated.

3.3

3.2.1 Metropolis method :

Assume that we have two states (o) and (n), and these states are
obeying Markov chain, where:
•
the N function denotes the probability density,
•
acc(initial state→final state) denotes the accepted move
from the initial state to final one,
•
푢
(Initial
state→final
state)
denotes
the transition
matrix
from the initial state to the final one and
푢 is the underlying
matrix of the Markov matrix, chosen to be symmetric.
In the equilibrium condition to move between the states is equal so
that the system does not change so:

3.4

The probability density of the certain state is the probability to find
this state inside the system,

3.5

where Z is configurational integral.

3.6

The choice of Metropolis :

3.8

3.2.2
Cluster moves technique :

The cluster moves technique is a sampling technique in order
to speed up the simulation. It is an artificial move and instead
of moving one particle from the cluster in an unfavorable
environment, one moves the particle with surrounding ions and
thereby one can speed up the simulation and so it increases the
acceptance rate. From the previous explanation of Metropolis
method, the acceptance of moving the cluster from the position (o)
to a new position (n),

But this equation doesn’t care about the equilibrium conditions
for the micelle with the counterions (the cluster in our model). In
order to modify the equation to serve the equilibrium conditions of
the cluster, the equation should not accept moving more than one
particle to or from the cluster. This term has been introduced p(k,l),
where p is the probability of moving the particle, which is a part
from the cluster between two positions k, l : p(k,l)= 1 when rk,l<
rc and = 0 when rk,l> rc
where rc is the radius of the cluster. The acceptance equation will
be

This equation guarantees the acceptance moves will obey the
equilibrium of the cluster.

3.2.3 Ensembles :

One special class of ensemble is the ones that do not evolve over
time. These ensembles are known as equilibrium ensembles and
their condition is known as statistical equilibrium. Statistical
equilibrium occurs if, for each state in the ensemble, the ensemble
also contains all of its future and past states with probabilities
equal to that state. The study of equilibrium ensembles of isolated
systems is the focus of statistical thermodynamics. In this study
the NVT ensemble has been used i.e. the number of molecules,
volume and temperature are being constant.

3.3
Structural analysis

A protein is a polymer formed naturally from a specific sequence of
amino acids (monomers), and the amino acids are bonded together
by peptide bonds, as shown in Figure 11. To calculate the size
of the proteins precisely, the radius of gyration is recommended
because it considers the distance between all the monomers, not
the average distance as in the hydrodynamic radius. The radius
will be discussed later in this section.

3.3.1
The radial distribution function

The radial distribution function (r.d.f or g(r)) is the variation of the number
density (𝜌 = !, the unit of 𝜌 is number/volume) of particle as a function of the
!
distance from a certain particle, see schematic picture in Figure 12 (left). Or more simply spoken: it is the probability of finding a particle at a distance r far from a reference particle.

The g(r) is calculated by using a Histogram technique, and it is a graphical representation of data distribution to estimate the probability distribution of a contentious variable. The radial distribution function is determined by calculating the distances between all the particles pairs, and put them all together in a histogram.

Figure 12: shows the radial distribution function g(r).[42] It is
a schematic picture to the left and the histogram in a graphi
-
cal representation to the right.

The graph of the radial distribution function can be understood
from the following: The beginning of moving from zero in g(r)
axis is the start of appearing a particles after the center particle,
and it is equal to twice of the radius, if it is a pure hard sphere
potential.
By assuming that (rmax) is the distance (r) at the maximum g(r),
means that at the (rmax) distance from the center particle is the
most crowded position. Also the height of the peak refers to the
number of the configuration. The g(r) gives a link between the
microscopic determined from the intermolecular forces, like in
this thesis ̈electrostatic interactions ̈, and the other microscopic
thermodynamic properties like energy (E) and entropy ( S).

3.3.2 End to end distance of the polymeric chain :

The distance from one end to another. By assuming that the
polymeric chain has N monomers, and bonds with length req
bond the polymers to each other. The end-
‐to- ‐end distance Ree
corresponds to the displacement Ld of a random walk with
number of steps Nr – 1 and step length req. Then

3.3.3
The hydrodynamic radius :

It is also related to the number of monomers and the distance
between monomers.

where N is number of monomers and rij the distance between i and
j monomers.

3.3.4The radius of gyration :

Is the root mean square of the distances between the monomers
and their centers of mass, and it can be expressed by the formula:

the ratio between the mean square of the end to end distance(Ree),
and the radius of gyration (Rg), is related to the factor, and the
shape factor is related to the polymer conformation.

4.Results and Discussion :

The charges of the micelles have been calculated according to
equation (2.22) at the pH of the normal fresh milk (6.7) and it
has been found to be (-
‐140e), which for 21 micelles gives rise
to 2940 counterions. In the simulations, all the counterions
have been explicitly taken into account, hence the systems are
electroneutral. The existence of charged species in the solution
gives rise to Columbic forces, i.e. attractive electrostatic
interactions between particles of the opposite sign, and repulsive
electrostatic interactions between particles of the same sign of
charge. Electrostatic interactions have an effect on the distribution
of micelles, and since the micelles repel each other, and attract the
counter ions, a cluster move technique has been applied i.e. it has
been possible to move a micelles with condensed counterions in
one MC step.

4.1The effect of the electrostatic forces :

The radial distribution function between the micelles is shown
in Figure 13 (blue curve), and the corresponding g(r) for non-
‐
charged micelles are shown in the red curve. The blue curve shows
the effect of the electrostatic interactions on the micellar structure.
The g(r) is zero for r less than one molecular diameter due to
hard sphere repulsive forces. The electrostatic interactions induce
repulsion between the micelles and there is a prominent peak at
around r = 328 Å, which is more than the double of the radius of
the micelles (150 Å). The height of the peak reaches two, which
means that it is two times more likely to find two micelles of this
distance than in an ideal gas.

The red curve shows the “ideal system ̈ with neutral micelles, and
as expected the micelles separating distance is the double of the
radius due to hard sphere interactions. In the ideal system g(r) = 1,
indicating no long-
‐ranged order.

By assuming a lattice model, and dividing the box volume by the
number of micelles, one can attain the maximum distance that can
be obtained between two micelles. In fact, this procedure gives a
distance of approximately 329 Å, hence the peak in Figure 13 (r
≈ 350 Å) seems to be a little bit off. This is due to a system size
artefact and Figure 14, shows the g(r) (blue curve) when the system
is made twice as big, i.e., the box volume and number of micellles
are doubled while keeping a constant micelle number density.
By increasing the volume, the correct inter-
‐particle micellar
distance is obtained. Due to the fact that the simulation time is
roughly increasing by a factor of N2 (N is the number of particles),
we have decided to use the smaller system in our investigation.
The obtained structural effects and qualitative trends are expected
to be similar.

4.2
The ideal system and the effect of counterion
volume

In Figure 15, the blue curve shows the “ideal system ̈ with
neutral micelles, and as expected, the closest approach distance
between the micelles is twice the micelle radius due to hard sphere
interactions. The red curve gives the structure for a system with
again 21 neutral micelles, but now including 2940 neutral hard
sphere “ions”. From this simple test, we conclude that the excluded
volume effect from the small ionic species is negligible.

Figure 15: The radial distribution function versus the distance
between the micelles for two different systems. The blue curve
corresponds to a system with 21 neutral micelles, while the red
curve corresponds to a system with 21 neutral micelles together
with 2940 neutral “counterions”.

4.3
The effect of salt concentration :

The monovalent salt concentration of 80 mM is interesting, since
it corresponds to the real ionic strength in milk. The salt solution
screens the electrostatic interactions and to be able to study the
effect of salt, a screened electrostatic potential has been used.
The g(r) in Figure 16 shows that when the salt concentration is
increased, the electrostatic interaction is not as strong as in the
reference system, and that it becomes more like the ideal system
c.f. green curve (ideal system) with the red curve (80 mM).

Figure 16: shows the response in the micelle-
‐micelle struc
-
ture with respect to the addition of 80 mM of monovalent
salt

The blue curve is for the reference system (explained above) and it
is clear that the micelles are more separated than in the red curve,
which is for the system with 80mM concentration (more screened
electrostatic interaction). The third system is the ideal solution
without any charges, where no long-
‐ranged order is visible. Notice
that the red and green curve deviate from each other at shorter
interparticle distances, hence at these distances the micelles “feel”
each other and repulsion is obtained. The screening length for 80
mM salt is ≈ 11 Å.

4.4
The effect of salt valence on milk :

The monovalent counterions in the reference system have been
replaced by divalent and tetravalent ions where the aim was to
study how the valency of the ions affects the structure, and if
correlation effects can be obtained. Correlation effects means that
two highly charged objects of the same sign can attract each other
due to charge fluctuations, see schematic mechanism in Figure
17. This contra intuitively effect is only visible for higher ordered
salts and most probably important in the aggregation process of
proteins

Figure 17: Correlations effects. Each half of the system is
electroneutral, and the blue spheres correspond to salt ions of
higher valency, and the surfaces correspond to negatively charged
colloids. In (a) the two colloids are well separated and do not feel
each other, (b) they feel each other and charge fluctuations are
induced, which results (c) in the aggregation of the colloids.

The blue curve in Figure 18 corresponds to the reference system
with monovalent counterions, the red curve corresponds to a
system where the monovalent ions have been replaced by divalent
ions, and the green curve corresponds to a system where the
monovalent ions have been replaced by tetravalent ions. As clearly
shown, when the valency is increasing the interparticle distance
between the micelles is decreasing and the largest effect is visible
for system with tetravalent ions. Here the micelles are more or
less in molecular contact with each other, only separated by the
counterions. Moreover, in the tetravalent system g(r) ≈ 13 at the
peak position, which indicates that at r
≈ 160 Å, it is 13 times higher probability to find neighboring
micelles in comparison with an ideal gas.

Figure 18: shows the difference in the distribution of micelles
with the variation of the valency of the counterions. In the
green curve, the micelles look freer than in the red and the blue
curve, and this behavior occurred because of the screening of
the electrostatic force, which increases sharply when increasing
the valence. The maxima of the r.d.f:s are varying as well,
which indicate the variations of structures.

4.5
The effect of pasteurization :

Increasing the temperature is of special interest, because that is
used in the pasteurization processes. Even without adding the
values of the right permittivity coefficient, the graph still shows
some variations of the behaviors of the micelles. These variations
occur because of the effect of temperature on the Debye screen
-
ing length, which is directly proportional to the root mean square
of the 푢!!.

Figure 20: The micellar structure has been
studied as an effect of temperature to mimic
the pasteurization process.

The radial distribution functions for the three different temperatures
are given in Figure 20. As shown, they coincide and one reason
might be because we kept the permittivity coefficient constant in
the simulation. At 740C, the permittivity coefficient of the water
falls to 63:61, and at 1380C it falls to 48:46.[40] The entropy
increases by increasing the temperature, but since the inter-
‐
micellar electrostatic repulsion is very strong, it does not have any
influence on the structure.

4.6
The effect of fermentation :

It is well-
‐known that the charge of the protein depends on the pH
as explained above. The protein charges have been calculated by
the equation (1.15) for different pH as shown in Figure 21

Figure 22: shows the g(r) for micelles in ordinary milk
(blue curve) and in fermented milk (red curve) as well as test micelle where the charge is 1/3 of the micellar charge
in milk.

As visible there are no traceable effects on micellar structure
between milk and fermented milk system (in which we increased
the micelle charge to 420 e.) This is due to the fact that the
screening of the counter decreases the Debye length to 21.5Å from
37.3Å for the reference system.

5.Summary :

Milk is indeed a complex liquid, which contains many different
components. In this study, we have used a simple model of
hard spheres to mimic
푢- ‐casein micelles in milk. The aim
was to investigate how electrostatic interactions affect the micellar
structure by varying the solution conditions.

The structure of the solution has been analyzed by comparing the
radial distribution function as a function of pH, salt concentration
and valency, and temperature. It was noticed that due to the fact
that the micellar charge is very large, the electrostatic repulsive
interaction dominates, and the mean distance between the micelles
are almost always obtained. The mean distance was calculated
to be approximately 330 Å for a micellar volume fraction of five
percent.

Moreover, an increase of the temperature does not affect the
structure at all i.e. the entropic contribution due to increased
temperature can be neglected in comparison with electrostatic
repulsion between the micelles.
When the salt concentration was increased to 80 mM, which
corresponds to the ionic strength in milk, the structure of the
푢- ‐casein micelles resembles the structure of an ideal gas i.e.
the electrostatic repulsive interactions are screened. 80 mM salt
corresponds to a screening length of approximately 11 Å

6.Suggestions for future work :

I.
Verifying the model by simulation larger systems.
II.
Investigate how a temperature affects the dielectric
constant, in our study; it was kept constant to 78.4.
III.
Investigated the effect of mixed micelles.
IV.
Investigate how the protein volume fraction affects the
structure.
V.
It is of interest to get a general understanding of the
micellar formation, and which theory that is applicable.

7.Acknowledgement :

I would like to express my deep gratitude to my master thesis
supervisors, Torbjörn Åkesson, and Marie Skepö, you have been a
tremendous mentor for me. I would like to thank you for letting a
part of this interesting project, and for your endless caring; Torbjörn
Åkesson: I will always remember you as a kind person, and I will
never forget our interesting discussions; Marie Skepö: Thank you
for giving me your time and care. Without Marie Skepö, this thesis
would not be that good.

Martin Turesson: Thank you for spending time to revise my thesis.
Henrik Dieden: from Klostergården Rotary Club, thank you for
your help and support.
Thank you all my previous professors.

Thank you all my colleges in theoretical chemistry department for
their support and friendship.
Thank you to my corridor mates, and my friends out side the
chemical center. Big thank you for my family for the wordless
love and support.