We wrote the perl script to calculate the center of mass. the equations are x=(m1*x1 + m2*x2...)/(m1+m2...), y=(m1*y1 + m2*y2...)/(m1+m2...) and z=(m1*z1 + m2*z2...)/(m1+m2...)
Two problems came up:
\$atoms -> CopyFrom(\$doc->DisplayRange->Sets(\$targetset)->Atoms); ------------------->we checked the xsd file and found it is empty.
for (my \$i=0; \$i< \$atoms->Atoms->Count; ++\$i) { ------------------------->infinite value/system. "Count" did not work. This problem might be triggered by the above one.
Please help~ here is the details of the perl script file:
--------------------------------------------------------------------------------------------------------------
#!perl
#The system includes a water droplet on a thin film in vacuum.
#This script is to calculate the center of Mass with a specific set (water).
use strict;
use MaterialsScript qw(:all);
###########################################################################
#Begin user editable settings
#
my \$filename = "filename"; #Name of the trajectory
my \$targetset = "water"; #Name of the target set
my \$surfaceposition = 32.25648; #Position of the surface
#
#End user editable settings
###########################################################################
#Defines the trajectory document
my \$doc = \$Documents{"\$filename.xtd"};
#create a new study table to hold the stractures and energies and set the headings
my \$newStudyTable = Documents->New("\$filename"."_CenterofMass.std");
my \$calcSheet = \$newStudyTable->Sheets->Item(0);
\$calcSheet->ColumnHeading(0) = "Frame Numer";
\$calcSheet->ColumnHeading(1) = "X";
\$calcSheet->ColumnHeading(2) = "Y";
\$calcSheet->ColumnHeading(3) = "Z";
#Get the total number of frames in the trajectory
my \$numFrames = \$doc->Trajectory->NumFrames;
print "Calculating the mass center for \$numFrames frames on \$filename trajectory\n";
#Initialise Forcite. If you want to change the settings, you can load them or
#use a change settings command here.
my \$forcite = Modules->Forcite;
#Starts to loop over the frames in the trajectory
for ( my \$count=1; \$count<=\$numFrames; \$count++) {
print "Calculating the center of mass for the frame \$count\n";
\$calcSheet->Cell(\$count-1, 0) = \$count;
#Define the frame of the trajectory
\$doc->Trajectorys->Item(0)->CurrentFrame = \$count;
#Create the target document with the target atoms
my \$atoms = Documents-> New("atoms.xsd");
\$atoms -> CopyFrom(\$doc->DisplayRange->Sets(\$targetset)->Atoms);
#Calculate the center of mass in the target atoms
my \$xcm = 0;
my \$ycm = 0;
my \$zcm = 0;
my \$mcm = 0;
for (my \$i=0; \$i< \$atoms->Atoms->Count; ++\$i) {
my \$atom = \$atoms->Atoms(\$i);
my \$myx = \$atom-> X;
my \$myy = \$atom-> Y;
my \$myz = \$atom-> Z;
my \$mymass = \$atom-> Mass;
if (\$myz >= \$surfaceposition) {
my \$Xmass = \$mymass*\$myx;
my \$Ymass = \$mymass*\$myy;
my \$Zmass = \$mymass*\$myz;
\$mcm = \$mcm + \$mymass;
\$xcm = \$xcm + \$Xmass;
\$ycm = \$ycm + \$Ymass;
\$zcm = \$zcm + \$Zmass;
\$Xmass->Discard;
\$Ymass->Discard;
\$Zmass->Discard;
}
else {
\$mcm = \$mcm;
\$xcm = \$xcm;
\$ycm = \$ycm;
\$zcm = \$zcm;
}
\$mymass->Discard;
\$myx->Discard;
\$myy->Discard;
\$myz->Discard;
\$atom->Discard;
}
my \$XCM = \$xcm/\$mcm;
my \$YCM = \$ycm/\$mcm;
my \$ZCM = \$zcm/\$mcm;
\$calcSheet->Cell(\$count-1, 1) = \$XCM;
\$calcSheet->Cell(\$count-1, 2) = \$YCM;
\$calcSheet->Cell(\$count-1, 3) = \$ZCM;
\$atoms->Discard;
\$XCM->Discard;
\$YCM->Discard;
\$ZCM->Discard;
}
print "Calculation complete.\n";