Commit 8e2c7ecd authored by Дарья Саликова's avatar Дарья Саликова
Browse files

Update md_statistics.jl

parent ed0d4d65
No related merge requests found
Showing with 17 additions and 5 deletions
+17 -5
......@@ -158,7 +158,7 @@ function c_mass(system::State{T}) where T
cmass .= cmass ./ molmass
return cmass
end
#автокоррелятор
function acf(input ::Vector{T}, acf_length ::Int64) where T
autocorr = zeros(T, acf_length)
for k in 1 : acf_length
......@@ -227,11 +227,12 @@ function msdt(direct::AbstractString, mask::AbstractString, dt::Integer)
systemref = MD.State()
molsyst = MD.State()
molref = MD.State()
print("------")
newmask = split(mask, '*')
newdir = readdir(direct)
accum = nothing
n_dumps = 0
for fname in readdir(direct) #выдаёт массив имён файлов, отсортировать массив строк
if startswith(fname, newmask[1]) && endswith(fname, newmask[2])
num = fname[sizeof(newmask[1])+1:end-sizeof(newmask[2])] #проверить, что они существуют прежде чем считывать
......@@ -261,12 +262,13 @@ function msdt(direct::AbstractString, mask::AbstractString, dt::Integer)
return accum / n_dumps
end
function rdfw!(gr::Vector{<:Real}, system::State, rmax::Real, nslice::Integer, weights::Vector{<:Real)
function rdfw!(gr::Vector{<:Real}, system::State, rmax::Real, nslice::Integer, weights::Vector{<:Real})
build_cell_list!(system, rmax)
coord = system.coord
cell_list = system.cell_list
dr = rmax / nslice
boxsize = system.size
ntotal = 0
V = prod(boxsize)
length(gr) == nslice || resize!(gr, nslice)
gr .= 0
......@@ -274,6 +276,7 @@ function rdfw!(gr::Vector{<:Real}, system::State, rmax::Real, nslice::Integer, w
neigh_i = neigh_atoms(system, i)
for atom1 in cell_list[i]
r1 = coord[atom1]
ntotal += weights[system.type[atom1]]
for atom2 in Iterators.reverse(neigh_i)
atom2 > atom1 || break
r2 = coord[atom2]
......@@ -285,7 +288,6 @@ function rdfw!(gr::Vector{<:Real}, system::State, rmax::Real, nslice::Integer, w
end
end
end
ntotal = natoms(system)
density = ntotal / V
for i in eachindex(gr)
vol_i = 4/3 * π * dr^3 * (i^3 - (i-1)^3)
......@@ -300,5 +302,15 @@ end
Compute radial distribution function with weights.
"""
rdfw(system::State, rmax::Real, nslice::Integer, weights::Vector{<:Real) = rdfw!(zeros(Float64, nslice), system, rmax, nslice, weights)
rdfw(system::State, rmax::Real, nslice::Integer, weights::Vector{<:Real}) = rdfw!(zeros(Float64, nslice), system, rmax, nslice, weights)
function autokorr_ave(input)
summ1 = 0
summ2 = 0
for idx in 1:2000000
a = readline(input)
summ1[idx] += tryparse.(Float64, a.split()[0])
summ2[idx] += tryparse.(Float64, a.split()[1])
end
print(summ1/2000000, summ2/2000000)
end
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment