Monte-Carlo in NVT ensemble for Lennard-Jonse potantioal in 2D
calculate the energy after 'Nstep' Monte-Carlo runs (Metropolis algorithm)
for a system of 'N' particles in 2D interacting via Lennard-Jonse pair
potential. periodic boundary conditions (PBC) are apllied.
inputs:
N - number of particles
T - reduced temperature
Nsteps - number of steps
maxdr - maximum particle displacement
initialConfig - initial configuration of particles (2 by N matrix)
initialU - initial energy of the configuration
rCutoff - the cutoff distance for the energy
optional: the 'm' in the interaction (U = 4((1/r)^12 - (1/r)^-m)),
default is 6.
outputs:
finalU - the energy in each step (1 by Nstep matrix)
finalConfigurations - the coordinates of all particles in each step (2 by N
by Nstep matrix)
finalDistances - the pair distances of all particles in each step (N by N
by Nstep matrix)
moveCount - counts accepted moves
the potential of Lennard-Jonse in reduced units:
U = 4*[(1/r)^12 - (1/r)^6]
(you can change the 6 with the optional input 'm')
the potential of Lennard-Jonse in non-reduced units:
U = 4*epsilon*[(sigma/r)^12 - (sigma/r)^6]
reduced units:
T(reduced) = kT/epsilon | r(reduced) = r/sigma | U(reduced) = U/epsilon
usage examples:
create an initial configuration of 100 particles in a 100*100 box:
N = 100;
L = 100;
initialConfig = L*rand(2,N);
rho = N/L^2;
choose simulation parameters:
T = 1; Nsteps = 1000; maxdr = 1; rCutoff = 2.5;
calculate the pair distances:
initialDistances = sqrt(bsxfun(@(x1,x2) (x1-x2).^2 ,...
initialConfig(1,:),initialConfig(1,:)')...
+bsxfun(@(x3,x4) (x4-x3).^2 ,...
initialConfig(2,:),initialConfig(2,:)'));
calculate the initial energy (up to the cutoff):
d = initialDistances(and(initialDistances < rCutoff,initialDistances > 0));
initialU = 4*sum(sum(d.^(-12)-d.^(-6)));
calculate energy
[finalU,~,~,finalConfiguration,finalDistances,moveCount] = ...
MonteCarlo2DLJHeart(N,T,rho,Nsteps,maxdr,initialConfig,rCutoff...
,initialDistances,initialU)
also calculate the pressure.
first calculate the initial virial:
virial = -(2*rho/N)*(sum(sum(6*d.^(-6) - 12*d.^(-12))));
now calculate the energy and pressure using monte carlo
[finalU,finalVirial,finalPressure,finalConfiguration,finalDistances,moveCount] = ...
MonteCarlo2DLJHeart(N,T,rho,Nsteps,maxdr,initialConfig,rCutoff...
,initialDistances,initialU,'virial',virial);
인용 양식
adi (2024). Monte Carlo simulation for two dimensional particles in the NVT ensemble-Lennard Jones interaction (https://www.mathworks.com/matlabcentral/fileexchange/55266-monte-carlo-simulation-for-two-dimensional-particles-in-the-nvt-ensemble-lennard-jones-interaction), MATLAB Central File Exchange. 검색됨 .
MATLAB 릴리스 호환 정보
플랫폼 호환성
Windows macOS Linux태그
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!버전 | 게시됨 | 릴리스 정보 | |
---|---|---|---|
1.0.0.0 |