Forum Navigation
You need to log in to create posts and topics.

MSD

Hello

i am trying to calculate the msd for particle #11 and tried this code:

import ovito
from ovito.io import import_file, export_file
from ovito.modifiers import *
import numpy

# Load input data and create a data pipeline
pipeline = import_file("Documents/results/dump.*")

# Calculate per-particle displacements with respect to initial simulation frame:
pipeline.modifiers.append(CalculateDisplacementsModifier())

# Define the custom modifier function:
def calculate_msd(frame, data):

# Access the per-particle displacement magnitudes computed by the
# CalculateDisplacementsModifier that precedes this user-defined modifier in the
# data pipeline:
displacement_magnitudes = data.particles['Displacement Magnitude']
particle_types = data.particles['Particle Type']
# Compute MSD:
for particle_index in range(data.particles.count):
if particle_types[particle_index] == 11:
msd = numpy.sum(displacement_magnitudes ** 2) / len(displacement_magnitudes)

# Output MSD value as a global attribute:
data.attributes["MSD"] = msd
data.attributes["Type"] = particle_types

# Insert user-defined modifier function into the data pipeline.
pipeline.modifiers.append(calculate_msd)

# Export calculated MSD value to a text file and let OVITO's data pipeline do the rest:
export_file(pipeline, "Documents/msd_data.txt",
format = "txt/attr",
columns = ["Timestep", "MSD"],
multiple_frames = True)

It seems to work but my MSD jumps in the first step from 0 to 0.4 and than keeps increasing linearly which i would expect from the beginning.

Hi Betim,

for the future, could you please make sure to use the "Code insert function" next to the "Insert/edit link" option to correctly format your python scripts when you post them here.

Please note that this part of your script probably does not do what you intended to do:

# Compute MSD:
for particle_index in range(data.particles.count):
    if particle_types[particle_index] == 11:
    msd = numpy.sum(displacement_magnitudes ** 2) / len(displacement_magnitudes)

Here, you loop over all particles, and then compute the total MSD (including all particle displacements) if the current particle is of type 11, which is a redundant operation.

In case you meant to calculate the MSD of those specific particles only, you can e.g. change your script to:

# Compute MSD and save as global attribute

data.attributes["MSD"] = numpy.sum(displacement_magnitudes ** 2) / len(displacement_magnitudes)
data.attributes["MSD Type 11"] = numpy.mean(displacement_magnitudes[particle_types == 11]**2)

 

Quote from Constanze Kalcher on August 25, 2020, 2:02 pm

Hi Betim,

for the future, could you please make sure to use the "Code insert function" next to the "Insert/edit link" option to correctly format your python scripts when you post them here.

Please note that this part of your script probably does not do what you intended to do:

# Compute MSD:
for particle_index in range(data.particles.count):
    if particle_types[particle_index] == 11:
    msd = numpy.sum(displacement_magnitudes ** 2) / len(displacement_magnitudes)

Here, you loop over all particles, and then compute the total MSD (inlcuding all particle displacements) if the current particle is of type 11, which is a redundant operation.

In case you meant to calculate the MSD of those specific particles only, you can e.g. change your script to:

# Compute MSD: 

data.attributes["MSD"] = numpy.sum(displacement_magnitudes ** 2) / len(displacement_magnitudes)
data.attributes["MSD Type 11"] = numpy.mean(displacement_magnitudes[particle_types == 11]**2)

 

Thanks a lot.