Add radial distribution plots and tests
Radial Distribution Functions
As previously, we calculate molecular RDFs for each distinct molecule type. Now, we also consider 4 distinct ranges of the trajectory for computing the ensemble average: first 20% of trajectory, last 80% of trajectory, last 60% of trajectory, last 40% of trajectory. In this way, one can see the convergence of the RDF throughout the MD workflow.
So, the GUI overview page should have a "Structural Properties" section, under which the RDFs are displayed. They should be organized by pairs of molecule types (following the visualization organization). I.e., each molecule type will have its own plot:
In terms of the location of quantities within the metainfo: archive.workflow[0].molecular_dynamics.results.radial_distribution_functions[0] holds all RDF info of a given "type". For us "type"==molecular and all the relevant quantities are in this section. There is a subsection called "radial_distribution_function_values" which is a list of the individual RDFs. "label" denotes the molecule-types involved in each RDF, "frame_start" and "frame_end" denote the start and stop time steps for the averaging, respectively, "bins" and "values" denote the x and y axes of the above plots, respectively.
Here is a code for grabbing the appropriate metainfo quantities and making the above plots:
section_MD = a.workflow[0].molecular_dynamics.results
sec_rdfs = section_MD.radial_distribution_functions[0]
pairs_list = []
for i_rdf, subsec_rdf in enumerate(sec_rdfs.radial_distribution_function_values):
pairs_list.append(subsec_rdf.label)
pairs_list = np.array(pairs_list)
pair_types = np.unique(pairs_list)
for i_pair_type, pair_type in enumerate(pair_types):
fig = plt.figure(figsize=(15,8))
indices = np.where(pairs_list == pair_type)[0]
for index in indices:
subsec_rdf = sec_rdfs.radial_distribution_function_values[index]
bins = ureg.convert(subsec_rdf.bins.magnitude, subsec_rdf.bins.units, ureg.angstrom)
values = subsec_rdf.value
frame_start = subsec_rdf.frame_start
frame_end = subsec_rdf.frame_end
plt.plot(bins, values, label='frames: '+str(frame_start)+'-'+str(frame_end), linewidth=3)
index = indices[-1]
subsec_rdf = sec_rdfs.radial_distribution_function_values[index]
plt.xlabel(r'r [$\AA$]', fontsize=22, labelpad=10)
plt.ylabel(r' $g(r)$', fontsize=22, labelpad=10)
if i_pair_type == 0:
plt.title(sec_rdfs.type+' radial distribution functions \n\n'+subsec_rdf.label, fontsize=22)
else:
plt.title(subsec_rdf.label, fontsize=22)
plt.xticks(fontsize=16)
plt.yticks(fontsize=16)
plt.legend(loc="best", fontsize=22)
plt.show()
Here are the related tasks:
-
New PropertyCard
for structural properties -
New component for RDF plots -
GUI tests -
Normalization tests for RDF