+1 (203) 503-8995


The header can only be two lines. The first one is not used but there must be some characters on it.

The second line states the file that you are going to be using (either .xls(x) or .pdb), the save file where the output comes and the gas that is being used. They must be separated by spaces (NO TABS).

Among  the  gas possibilities you have the following:

N2, He, Argon, CO2, SF6, Air

The code is case sensitive in this case.The selection of the gas here only affects which LJtable.xlsx sheet is used for the Lennard-Jones parameters.



File                                      Savefile                  Gas

pdb1y32mod12.pdb       Saveme.txt              N2


The body consists of the parameters that can be altered in the code. Although some of these parameters can be omitted in this file, I would suggest keeping them all for future references. We will describe parameter by parameter. Note that the parameter and its value are separated by a space.


The first parameter is [interface]. It has two numbers separated by spaces. Both values must be there for the program to work although the second value is not used at the moment. The first value states (Windows only for now) whether you want the interface to show up or not. Setting a value of 1 will open the interface while setting a value of 0 will make the command window popup but not the interface. We will call this last mode, GUIFree.

interface 0 0

The rest of the parameters herein are only used in the GUIFree mode.

[tovalue] and [fromvalue] are used to establish the excelsheets that you want the program to compute in your excelfile. If a .pdb file is chosen, it will calculate more than one model from the .pdb file if more than one model exists. If there is not more than one model, make sure that both values are equal so that the program does not repeat the same calculation several times. In the .pdb file case, it does not matter what [fromvalue] is chosen, it will always start from model 1. The way to establish several models in a .pdb file is as follows:

  • You must have a TER line to terminate a chain in your model
  • You must have an ENDMDL line to terminate your model
  • You must have a MODEL line that starts your new model

Example from Humanin (file : pdb1y32mod12.pdb)

ATOM    390  HB1 ALA A  24       1.791  10.249 -16.257  1.00  0.00           H

ATOM    391  HB2 ALA A  24       0.371   9.436 -16.916  1.00  0.00           H

ATOM    392  HB3 ALA A  24       1.609   8.504 -16.072  1.00  0.00           H

TER     393      ALA A  24


MODEL       13

ATOM      1  N   MET A   1       8.093  -8.495  15.332  1.00  0.00           N

ATOM      2  CA  MET A   1       6.757  -8.803  15.906  1.00  0.00           C

ATOM      3  C   MET A   1       5.730  -7.748  15.511  1.00  0.00           C

fromvalue 1

tovalue 2

The program takes the [charge] value set here in the following cases:

  • You have not assigned a total charge in your excel file regardless of how many charges have been assigned to atoms.
  • Any .pdb file regardless of how many charges have been assigned.

If no charge is assigned to any atoms, the charge established here will be set in the geometrical center of the structure. IT can be a positive or negative value.

Charge 1

[Mgas] value establishes the mass of gas used in Da. Make sure it agrees with your choice of gas on the header

Mgas 28

[radgas] is the value of the gas radii in Angstroms. Normally it is around 1.5 for N2 and 1.2 for He. These values are used for PA methods and any hard sphere methods.

radgas 1.5

[Polarizability] is the polarizability of the gas in A3. Normally 1.7A3 for N2 and 0.2073A3 for He.

Polarizability 1.7

[Pressure] is the pressure in Pascals.

Pressure 101325

[Mweight] is the molecular weight of your molecule in Da. It is only read here in the same cases as [charge]. For any .pdb file and in case it is not assigned in the excel file.

Mweight 256000

[Temperature] is the temperature used in K. This Temperature is only used for the speed of the gas and the calculation of the mobility. It does not affect the structure in any way.

Temperature 304

[redcoef] reduces the size of the  structure by the percentage applied in this coefficient. A value of 1.2 will reduce the dimensions of the distance to the geometric center and the radii of the atoms by 20%.

redCoef 1.0

[NrotationsPA] establishes the number of random orientations done in the Projected Area algorithm . Not used if [PA] is set to 0.

NrotationsPA 100

[NrotationsEHSS] and [NrotationsTM] establish the number of times the molecule is rotated when doing these calculations. If [Mobcal] is set to 0, these values must be set to 3. This will calculate the drag in 3 perpendicular directions and find the pral. Directions. If [Mobcal] is set to 1, it will operate in a similar fashion to Mobcal algorithms where it will pick up random directions and do the average value over all orientations. When [Mobcal]=1, normally [NrotationsEHSS]=[NrotationsTM]=30 is sufficient with [NgastotalEHSS] and [NgastotalTM] being  3000 and 5000 respectively.

NrotationsEHSS 3

NrotationsTM 3

[NgastotalEHSS] and [NgastotalTM] establish the number of gas molecules per rotation. If [Mobcal]=0, a value of 30000 for EHSS is normally enough and a value of 30000 to 50000 is normally enough for TM. The greater this value, the more accurate the calculation will be. If [Mobcal] =1 refer to the [NrotationsEHSS] explanation.

NgastotalEHSS 30000

NgastotalTM 30000

[Acommodation] is the accommodation coefficient for the EHSS/DHSS , TM and DTM methods that do not use Lennard-Jones (all those that have HS collisions). It has a value between 0 and 1 and describes which percentage of the collisions are purely specular and elastic and which are not. If 0 is chosen, all collisions will be specular and elastic (EHSS). If 1 is chosen, 100% of the collision will not, in principle, be specular and elastic.

The [diffuse?], [reemvel] and [other] boxes below only act upon the percentage that is not fully specular and elastic, so if [Acommodation]  is chosen as 0.3, 30% will be affected by these boxes. The [diffuse?] box tells you if the re-emission angles of that 30% are specular (set to 0) or diffuse (set to 1). The re-emission velocities [reemvel] will give you the inelasticity of the collision.

Acommodation 1.0

[Timestep] is only used for Trajectory methods and is used to increase or decrease the timestep in the velocity verlet algorithm. The higher this value, the smaller the timestep is. Values should not be changed unless tried before. Normally, heavily charged species require higher values of [Timestep].

Timestep 150

[Boxdomain] should not be changed

Boxdomain 16

[Diffuse?] has a value of 0 or 1 and it will tell you if the value given by [accommodation] will be reemitted as diffuse([Diffuse]=1) or only as specular([Diffuse]=0).

Diffuse? 1

The reemission velocities [reemvel] can be chosen to be taken from a Maxwell distribution (distributed Maxwell , distributed 92% Maxwell ,distributed Other ) or they can be taken as the mean values of such distributions which is recommended only for EHSS/DHSS (Maxwell, 92% Maxwell, Other). Choosing elastic means that you will re-emit with the same velocity as the gas molecule impinged and choosing elastic(mean) tells you that your gas molecules will be re-emitted with the mean of the incoming velocities taken from a skewed Maxwellian distribution. The [reemvel] is a numerical value. These are the corresponding velocities to the numbers used:

1 Elastic(mean), 2 Elastic, 3 Maxwell, 4 92%Maxwell, 5 Maxwell Distrib, 6 92% Maxwell Distrib, 7 Other, 8 Other Distrib

reemvel 4

If [reemvel] is given a value of 7 or 8, [other] is used to establish the mean value of the Maxwellian distribution employed.

Other 0

[Simplify] will simplify your structure so that you can get rid of all the internal atoms that will not affect Hard Sphere collision. It is very recommended that this value is set to 1 unless a [LennardJones] is set to 1 for TM methods.

Simplify 0

[PA],[EHSS/DHSS],[TM] and [DTM] have only values of 0 or 1 and establish which calculations will be performed.  You can perform more than one calculation at a time.

PA 1


TM 0


[PATSA] establishes which version of the PA algorithm is used. If [PATSA]=1, a Traveler Salesman(TSA) type of calculation is performed. IF [PATSA]=0 a MonteCarlo(MC) approach is used. TSA is about 30 times faster so it is recommended in general. If TSA causes problems or there are passable holes in the structure, the MC approach should be used instead.


[Cutoff] is only used if [LennardJones] methods are used and it establishes a cutoff value. It is not recommended unless the structure is very large (>10000 atoms).

Cutoff 0

[Mobcal] establishes whether or not principal directions are used. If [Mobcal] =0, then [NrotationsEHSS] and [NrotationsTM] must be set to 3 and pral. Directions will be calculated. If [Mobcal]=1, then a Mobcal approach where random directions are used and averaged. See [NrotationsEHSS] for further explanation.

Mobcal 0

[LennardJones] is only used in TM methods. If [LennardJones]=0, a 4-inf potential will be used when [TM] is selected. If [LennardJones]=1, LJ parameters  in LJTable will be used and a 4-6-12 potential will be used in [TM].

LennardJones 1

[SimplifiedTM] tells whether a simplification of the collision integral is used or not. In general 10 to 100 times more gas molecules must be used if [SimplifiedTM]=0 to achieve similar accuracy. It is therefore not recommended at the moment and should be set to [SimplifiedTM]=1. See Larriba and Hogan JCP 2013 for more information.

SimplifiedTM 1






Explanation of the EXCEL FILE

% Developed by Carlos LARRIBA ANDALUZ
% sometime in 2011-2013

% This routine forms part of a much larger project of calculating the %
% mobility of molecular clusters in a gas. In this case, the mobility is %
% calculated by approximating the momentum transfer of gas molecules as %
% their trajectory is bent by the action of an induced dipole. This %
% momentum transfer can be related to the equation for mobility: %
% Z= (1+zeta*Z^2)*3/16/Chi*z/p/Omega*(2*pi*kb*T/m_gas)^.5
% Here, %
% Z is the electrical mobility %
% z is the charge %
% p is the pressure in the DMA %
% Omega is the collisional cross section %
% kb is boltzmann’s constant %
% T is the temperature %
% m_gas is the molecular mass of the gas %

%%%%%%%%%%%%%%%%%%%%%%%%%%%%% EXCEL FILES %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Excel files should be formatted as follows:
% Choose an excelsheet and give it a name. Several excel sheets can be
% used on the same excel file.
% First column can be any string as long as it is not a number.
% In general, this column is reserved for the atom label. It is not
% necessary. The progam will obviate this column as long as there are
% no numbers on it.
% Columns 2-4 are the x, y, and z coordinates of the center of your atom.
% The units are Angstroms.
% Column 5 is the radii of the atoms in Angstroms.
% Column 6 is the charge of the atoms. Blank if not charged. The charges
% can be non-integers,i.e partial charges.
% Second line of Column 7 can be used to specify the total number of charges
% in your molecule.
% Column 8 specifies the Lennard-Jones parameters and the color of the ion depending on its amu(rounded to
% the nearest integer). Atoms with available Lennard-Jones parameters can be viewed in the file LJtable.xlsx

% Any data after Column 7 will be neglected.
% Sample file C60positions.xlsx is available for clarifications