diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index 2de06c5..fe46847 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -16,6 +16,6 @@ jobs: - name: Make sure virtualevn>20 is installed, which will yield newer pip and possibility to pin pip version. run: pip install "virtualenv>20" - name: Install Tox - run: pip install tox + run: pip install tox - name: Run pre-commit in Tox run: tox -e pre-commit diff --git a/aiida_uppasd/tools/__init__.py b/aiida_uppasd/tools/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/aiida_uppasd/tools/core_cli.py b/aiida_uppasd/tools/core_cli.py index 0635fc8..743d0d0 100644 --- a/aiida_uppasd/tools/core_cli.py +++ b/aiida_uppasd/tools/core_cli.py @@ -26,7 +26,7 @@ def uppasd_cli(): @click.option('-plot_name', default='None') @click.option('-width', default=100) @click.option('-height', default=20) -@click.argument('pk') +@click.argument('node_pk') def visualization_observations( iter_slice: int, y_axis: list, @@ -34,7 +34,7 @@ def visualization_observations( plot_name: str, width: int, height: int, - pk: int, + node_pk: int, ): # pylint: disable=too-many-arguments, too-many-branches, too-many-statements """ Visualize a given observable @@ -51,11 +51,11 @@ def visualization_observations( :type width: int :param height: height of the figure :type height: int - :param pk: pk number of the calculation that one wishes to visualize - :type pk: int + :param node_pk: node_pk number of the calculation that one wishes to visualize + :type node_pk: int """ - auto_name = locals() - cal_node = orm.load_node(pk) + auto_name = {} + cal_node = orm.load_node(node_pk) if iter_slice != -1: for name in y_axis: @@ -66,7 +66,7 @@ def visualization_observations( iter_list = cal_node.get_array('iterations')[:int(iter_slice)].astype(int) for name in y_axis: _name = str(name) - plotext.plot(iter_list, eval(_name), label=_name) + plotext.plot(iter_list, _name, label=_name) plotext.plotsize(width, height) if plot_name != 'None': plotext.title(f'{plot_name}') @@ -77,7 +77,7 @@ def visualization_observations( iter_list = cal_node.get_array('iterations')[:int(iter_slice)].astype(int) for name in y_axis: _name = str(name) - plotext.scatter(iter_list, eval(_name), label=_name) + plotext.scatter(iter_list, _name, label=_name) plotext.plotsize(width, height) if plot_name != 'None': plotext.title(f'{plot_name}') @@ -96,7 +96,7 @@ def visualization_observations( iter_list = cal_node.get_array('iterations').astype(int) for name in y_axis: _name = str(name) - plotext.plot(iter_list, eval(_name), label=_name) + plotext.plot(iter_list, _name, label=_name) plotext.plotsize(width, height) if plot_name != 'None': plotext.title(f'{plot_name}') @@ -107,7 +107,7 @@ def visualization_observations( iter_list = cal_node.get_array('iterations')[:int(iter_slice)].astype(int) for name in y_axis: _name = str(name) - plt.scatter(iter_list, eval(_name), label=_name) + plt.scatter(iter_list, _name, label=_name) plotext.plotsize(width, height) if plot_name != 'None': plotext.title(f'{plot_name}') @@ -119,15 +119,15 @@ def visualization_observations( def output_node_query( - cal_node_pk: typing.Union[int, str], + cal_node_node_pk: typing.Union[int, str], output_array_name: str, attribute_name: str, ) -> np.ndarray: """ Get the array output of a given calculation node - :param cal_node_pk: pk number for the node that is being queried - :type cal_node_pk: typing.Union[int, str] + :param cal_node_node_pk: node_pk number for the node that is being queried + :type cal_node_node_pk: typing.Union[int, str] :param output_array_name: name of the array that we are looking for :type output_array_name: str :param attribute_name: specific entry in the array that one is looking for @@ -138,7 +138,7 @@ def output_node_query( query_builder = orm.QueryBuilder() query_builder.append( orm.CalcJobNode, - filters={'id': str(cal_node_pk)}, + filters={'id': str(cal_node_node_pk)}, tag='cal_node', ) query_builder.append( @@ -178,7 +178,7 @@ def trajectory_parser( return mom_states -def get_arrow_next(array_name: str): +def get_arrow_next(array_name: str, rotation, mom_array_from_result, coord_r): """ Get the coordinates of an array in cartesian and in rotated coordinates. @@ -187,7 +187,7 @@ def get_arrow_next(array_name: str): :return: coordinates in cartesian and rotated coordinates :rtype: typing.Union[np.array,np.array,np.array,np.array, np.array, np.array] """ - rot_mom_array = r.apply(mom_array_from_result[array_name]) + rot_mom_array = rotation.apply(mom_array_from_result[array_name]) return ( coord_r[:, 0], coord_r[:, 1], @@ -198,14 +198,24 @@ def get_arrow_next(array_name: str): ) -def animate(data): +def animate( + quivers, + axes, + array_name, + rotation, + mom_array_from_result, + coord_r, + arrow_ratio_arr, + length_arr, + colors_arr, + normalize_flag_arr, +): #pylint: disable=too-many-arguments """ Animate the magnetic moments """ - global quivers quivers.remove() - quivers = ax.quiver( - *get_arrow_next(data), + quivers = axes.quiver( + *get_arrow_next(array_name, rotation, mom_array_from_result, coord_r), arrow_length_ratio=arrow_ratio_arr, length=length_arr, colors=colors_arr, @@ -222,14 +232,13 @@ def animate(data): @click.option('-normalize_flag', default=True) @click.option('-height', default=20) @click.option('-width', default=20) -@click.option('-color_bar_axis', default='x') @click.option('-path_animation', default='./motion.gif') @click.option('-interval_time', default=200) @click.option('-dpi_setting', default=100) @click.option('-path_frame', default='./motion.png') @click.option('-frame_number', default=0) @click.option('-animation_flag', default=False) -@click.argument('pk') +@click.argument('node_pk') def visualization_motion( rotation_axis: str, rotation_matrix: list, @@ -239,14 +248,13 @@ def visualization_motion( normalize_flag: bool, height: int, width: int, - color_bar_axis: str, path_animation: str, interval_time: int, dpi_setting: int, path_frame: str, frame_number: int, animation_flag: bool, - pk: int, + node_pk: int, ): # pylint: disable=too-many-arguments, too-many-locals """ Visualize the magnetic moments of the calculation node @@ -267,8 +275,6 @@ def visualization_motion( :type height: int :param width: width of the plot :type width: int - :param color_bar_axis: which axis determines the color bar - :type color_bar_axis: str :param path_animation: path to store the animation :type path_animation: str :param interval_time: how often one performs the animation @@ -281,18 +287,16 @@ def visualization_motion( :type frame_number: int :param animation_flag: whether or not to animate the figure :type animation_flag: bool - :param pk: pk number of the calculation that one wishes to animate the moments - :type pk: int + :param node_pk: node_pk number of the calculation that one wishes to animate the moments + :type node_pk: int """ - global coord_r, mom_array_from_result, r, quivers, axis_to_colorbar, ax - global arrow_ratio_arr, length_arr, colors_arr, normalize_flag_arr - r = Rotation.from_euler(rotation_axis, rotation_matrix, degrees=True) - mom_states_x = output_node_query(pk, 'trajectories_moments', 'moments_x') - mom_states_y = output_node_query(pk, 'trajectories_moments', 'moments_y') - mom_states_z = output_node_query(pk, 'trajectories_moments', 'moments_z') - coord = output_node_query(pk, 'coord', 'coord')[:, 1:4] + rotation = Rotation.from_euler(rotation_axis, rotation_matrix, degrees=True) + mom_states_x = output_node_query(node_pk, 'trajectories_moments', 'moments_x') + mom_states_y = output_node_query(node_pk, 'trajectories_moments', 'moments_y') + mom_states_z = output_node_query(node_pk, 'trajectories_moments', 'moments_z') + coord = output_node_query(node_pk, 'coord', 'coord')[:, 1:4] - coord_r = r.apply(coord) + coord_r = rotation.apply(coord) atoms_total = len(coord) arrow_ratio_arr = arrow_head_ratio length_arr = length_ratio @@ -305,19 +309,12 @@ def visualization_motion( atoms_total, ) - if color_bar_axis == 'x': - axis_to_colorbar = 0 - elif color_bar_axis == 'y': - axis_to_colorbar = 1 - else: - axis_to_colorbar = 2 - fig = plt.figure(figsize=(height, width)) - ax = fig.gca(projection='3d') + axes = fig.gca(projection='3d') if not animation_flag: - quivers = ax.quiver( - *get_arrow_next(frame_number), + quivers = axes.quiver( + *get_arrow_next(frame_number, rotation, mom_array_from_result, coord_r), arrow_length_ratio=arrow_ratio_arr, length=length_arr, colors=colors_arr, @@ -326,8 +323,8 @@ def visualization_motion( fig.savefig(path_frame) if animation_flag: - quivers = ax.quiver( - *get_arrow_next(0), + quivers = axes.quiver( + *get_arrow_next(0, rotation, mom_array_from_result, coord_r), arrow_length_ratio=arrow_ratio_arr, length=length_arr, colors=colors_arr, @@ -336,6 +333,10 @@ def visualization_motion( ani = FuncAnimation( fig, animate, + fargs=( + quivers, axes, rotation, mom_array_from_result, coord_r, arrow_ratio_arr, length_arr, colors_arr, + normalize_flag_arr + ), frames=list(range(len(mom_array_from_result))), interval=interval_time, ) diff --git a/aiida_uppasd/tools/post_processing/__init__.py b/aiida_uppasd/tools/post_processing/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/aiida_uppasd/tools/post_processing/plotting.py b/aiida_uppasd/tools/post_processing/plotting.py new file mode 100644 index 0000000..639ee8e --- /dev/null +++ b/aiida_uppasd/tools/post_processing/plotting.py @@ -0,0 +1,920 @@ +"""Set of functions to plot uppasd results""" +import os +from typing import Union + +from PIL import Image +import matplotlib.pyplot as plt +import numpy as np +import plotly.graph_objects as go +from scipy import ndimage +import vtk + +from aiida import orm + + +def plot_skyrmion_number(node: orm.Node, plot_dir: Union[str, os.PathLike], method: str = 'plain'): + """Plot the skyrmion number""" + + if method.lower() in ['plain', 'heatmap_plain']: + _axis = 0 + if method.lower() == ['average', 'heatmap_average']: + _axis = 1 + + if method.lower() not in ['plain', 'average', 'heatmap_plain', 'heatmap_average']: + raise ValueError('The method should be either "plain" or "average"') + + _sk_data = [] + _bfield = [] + _temperature = [] + for _sub_node in node.called_descendants: + _sk_data.append( + _sub_node.get_outgoing().get_node_by_label('sk_num_out').get_array('sk_number_data')[0][:, _axis] + ) + _bfield.append(_sub_node.inputs.inpsd_dict['hfield'].split()[-1]) + _temperature.append(_sub_node.inputs.inpsd_dict['temp']) + + if method.lower() in ['plain', 'average']: + _bfield_unique = np.unique(_bfield)[0] + + fig = go.Figure() + + for _field in _bfield_unique: + _temperature_curr = np.asarray(_temperature)[np.where(np.asarray(_bfield) == _field)[0]] + _sk_number_curr = np.asarray(_sk_data[np.where(np.asarray(_bfield) == _field)[0]]) + + fig.add_trace(go.Scatter(_temperature_curr, _sk_number_curr, name=f'B: {_field} T')) + fig.update_layout( + xaxis_title='Temperature', + yaxis_title='Skyrmion number', + ) + + fig.write_image(f'{plot_dir}/sk_number_{method}.png') + + if method.lower() in ['heatmap_plain', 'heatmap_average']: + + fig = go.Figure(data=go.Heatmap( + z=_sk_data, + x=_temperature, + y=_bfield, + zsmooth='best', + zmid=0, + )) + + fig.update_xaxes(range=[min(_temperature), max(_temperature)]) + fig.update_yaxes(range=[min(_bfield), max(_bfield)]) + fig.update_traces(colorscale='Jet', selector=dict(type='heatmap')) + fig.write_image(f'{plot_dir}/sk_number_heatmap.png') + + +def combined_for_pd(node: orm.Node, plot_dir: Union[str, os.PathLike]): + """Combine the collected data for the phase diagram""" + if node.inputs.plot_individual.value > int(0) and node.inputs.plot_combine.value > int(0): + skyrmions_singleplot = [] + for _index_temperature, _ in enumerate(node.inputs.temperatures): + for _index_field, _ in enumerate(node.inputs.external_fields): + skyrmions_singleplot.append( + Image.open(f'{plot_dir}/T{_index_temperature}B{_index_field}.png').crop((200, 0, 2096, 1280)) + ) #we need to crop the white edges + phase_diagram = Image.new('RGB', (648 * len(node.inputs.temperatures), 640 * len(node.inputs.external_fields))) + for _index_temperature, _ in enumerate(node.inputs.temperatures): + for _index_field, _ in enumerate(node.inputs.external_fields): + phase_diagram.paste( + skyrmions_singleplot[len(node.inputs.external_fields) * _index_temperature + _index_field], + (0 + 648 * _index_temperature, 0 + 640 * _index_field) + ) + phase_diagram.save(f'{plot_dir}/PhaseDiagram.png', quality=100) + + +# Read Location of Atoms +def read_atoms(coord_array_data, max_number_atoms): + """Transform the atomic positions to VTK structure""" + + points = vtk.vtkPoints() + number_atoms = 0 + # Read ahead + + # Read all data + for data in coord_array_data.get_array('coord'): + if number_atoms <= max_number_atoms: + coord_x, coord_y, coord_z = float(data[1]), float(data[2]), float(data[3]) + points.InsertPoint(number_atoms, coord_x, coord_y, coord_z) + number_atoms = number_atoms + 1 + return points, number_atoms + + +# Read vectors +# We must know the time step and the number of atoms per time +def read_vectors_data(mom_states, number_atoms): + """Transform the magnetic moments to VTK structure""" + final_mom_states_x = mom_states.get_array('mom_states_x')[-number_atoms:] + final_mom_states_y = mom_states.get_array('mom_states_y')[-number_atoms:] + final_mom_states_z = mom_states.get_array('mom_states_z')[-number_atoms:] + # Create a Double array which represents the vectors + vectors = vtk.vtkFloatArray() + colors = vtk.vtkFloatArray() + + # Define number of elemnts + vectors.SetNumberOfComponents(3) + colors.SetNumberOfComponents(1) + + for index in range(number_atoms): + vec_x, vec_y, vec_z = float(final_mom_states_x[index]), float(final_mom_states_y[index] + ), float(final_mom_states_z[index]) + _color = (vec_z + 1.0) / 2.0 + vectors.InsertTuple3(index, vec_x, vec_y, vec_z) + colors.InsertValue(index, _color) + return vectors, colors + + +#---------------------------- +# screenshot code begins here +#---------------------------- +#A function that takes a renderwindow and saves its contents to a .png file +def screenshot(ren_win, temperature, field, plot_dir): + """Take a screenshot of the magnetic configuration and save it as a png""" + win2im = vtk.vtkWindowToImageFilter() + win2im.ReadFrontBufferOff() + win2im.SetInput(ren_win) + # + povexp = vtk.vtkPOVExporter() + povexp.SetRenderWindow(ren_win) + #povexp.SetInput(renWin) + ren_win.Render() + # ren_win.SetFileName('realspace_plot_T_{}_B_{}.pov'.format(T,B)) + # povexp.Write() + # + to_png = vtk.vtkPNGWriter() + to_png.SetFileName(f'{plot_dir}/T{temperature}B{field}.png') + to_png.SetInputConnection(win2im.GetOutputPort()) + to_png.Write() + return 0 + + +def plot_realspace(points, vectors, colors, temperature, field, plot_dir): #pylint: disable=too-many-statements,too-many-locals + """Plot the real space magnetic configuration""" + ren_win = vtk.vtkRenderWindow() + ren = vtk.vtkOpenGLRenderer() + ren_win.AddRenderer(ren) + + # Set color of backgroung + ren.SetBackground(1.0, 1.0, 1.0) + ren_win.SetSize(2096, 1280) + + data_test = vtk.vtkPolyData() + data_scal = vtk.vtkUnstructuredGrid() + + # Read atom positions + atom_data = points + data_test.SetPoints(atom_data) + data_scal.SetPoints(atom_data) + + # Read data for vectors + vecz, colz = (vectors, colors) + data_test.GetPointData().SetVectors(vecz) + data_test.GetPointData().SetScalars(colz) + data_scal.GetPointData().SetScalars(colz) + + # Create colortable for the coloring of the vectors + lut = vtk.vtkLookupTable() + for i in range(0, 128, 1): + lut.SetTableValue(i, (127.0 - i) / 127.0, i / 127.0, 0, 1) + for i in range(128, 256, 1): + lut.SetTableValue(i, 0, (256.0 - i) / 128.0, (i - 128.0) / 128, 1) + lut.SetTableRange(-1.0, 1.0) + lut.Build() + + # Set up atoms + ball = vtk.vtkSphereSource() + ball.SetRadius(1.00) + ball.SetThetaResolution(16) + ball.SetPhiResolution(16) + + balls = vtk.vtkGlyph3DMapper() + balls.SetInputData(data_test) + balls.SetSourceConnection(ball.GetOutputPort()) + balls.SetScaleFactor(0.15) + balls.SetScaleModeToNodata_scaling() + balls.SetLookupTable(lut) + balls.Update() + + atom = vtk.vtkLODActor() + atom.SetMapper(balls) + atom.GetProperty().SetOpacity(1.5) + xmin, xmax = atom.GetXRange() + ymin, ymax = atom.GetYRange() + zmin, zmax = atom.GetZRange() + xmid = (xmin + xmax) / 2 + ymid = (ymin + ymax) / 2 + zmid = (zmin + zmax) / 2 + atom.SetPosition(-xmid, -ymid, -zmid) + atom.GetProperty().SetSpecular(0.3) + atom.GetProperty().SetSpecularPower(60) + atom.GetProperty().SetAmbient(0.2) + atom.GetProperty().SetDiffuse(0.8) + + # Create vectors + arrow = vtk.vtkArrowSource() + arrow.SetTipRadius(0.25) + arrow.SetShaftRadius(0.15) + arrow.SetTipResolution(72) + arrow.SetShaftResolution(72) + + glyph = vtk.vtkGlyph3DMapper() + glyph.SetSourceConnection(arrow.GetOutputPort()) + glyph.SetInputData(data_test) + glyph.SetScaleFactor(2.00) + # Color the vectors according to magnitude + glyph.SetScaleModeToNodata_scaling() + glyph.SetLookupTable(lut) + glyph.SetColorModeToMapScalars() + glyph.Update() + + vector = vtk.vtkLODActor() + vector.SetMapper(glyph) + vector.SetPosition(-xmid, -ymid, -zmid) + + vector.GetProperty().SetSpecular(0.3) + vector.GetProperty().SetSpecularPower(60) + vector.GetProperty().SetAmbient(0.2) + vector.GetProperty().SetDiffuse(0.8) + vector.GetProperty().SetOpacity(1.0) + + # Scalar bar + scalar_bar = vtk.vtkScalarBarActor() + scalar_bar.SetLookupTable(lut) + scalar_bar.SetOrientationToHorizontal() + scalar_bar.SetNumberOfLabels(0) + scalar_bar.SetPosition(0.1, 0.05) + scalar_bar.SetWidth(0.85) + scalar_bar.SetHeight(0.3) + scalar_bar.GetLabelTextProperty().SetFontSize(8) + + #Depth sorted field + data_test = vtk.vtkDepthSortPolyData() + data_test.SetInputData(data_test) + data_test.SetCamera(ren.GetActiveCamera()) + data_test.SortScalarsOn() + data_test.Update() + + #cubes + cube = vtk.vtkCubeSource() + cubes = vtk.vtkGlyph3DMapper() + cubes.SetInputConnection(data_test.GetOutputPort()) + cubes.SetSourceConnection(cube.GetOutputPort()) + cubes.SetScaleModeToNodata_scaling() + cubes.SetScaleFactor(0.995) + cubes.SetLookupTable(lut) + cubes.SetColorModeToMapScalars() + cubes.ScalarVisibilityOn() + cubes.OrientOff() + cubes.Update() + cubes.SetLookupTable(lut) + + cube_actor = vtk.vtkActor() + cube_actor.SetMapper(cubes) + cube_actor.GetProperty().SetOpacity(0.05) + cube_actor.SetPosition(-xmid, -ymid, -zmid) + + #hedgehog + hhog = vtk.vtkHedgeHog() + hhog.SetInputData(data_test) + hhog.SetScaleFactor(5.0) + hhog_mapper = vtk.vtkPolyDataMapper() + hhog_mapper.SetInputConnection(hhog.GetOutputPort()) + hhog_mapper.SetLookupTable(lut) + hhog_mapper.ScalarVisibilityOn() + hhog_actor = vtk.vtkActor() + hhog_actor.SetMapper(hhog_mapper) + hhog_actor.SetPosition(-xmid, -ymid, -zmid) + + #cut plane + plane = vtk.vtkPlane() + plane.SetOrigin(data_scal.GetCenter()) + plane.SetNormal(0.0, 0.0, 1.0) + plane_cut = vtk.vtkCutter() + plane_cut.SetInputData(data_scal) + plane_cut.SetInputArrayToProcess( + 0, 0, 0, + vtk.vtkDataObject().FIELD_ASSOCIATION_POINTS, + vtk.vtkDataSetAttributes().SCALARS + ) + plane_cut.SetCutFunction(plane) + plane_cut.GenerateCutScalarsOff() + plane_cut.SetSortByToSortByCell() + clut = vtk.vtkLookupTable() + clut.SetHueRange(0, .67) + clut.Build() + plane_mapper = vtk.vtkPolyDataMapper() + plane_mapper.SetInputConnection(plane_cut.GetOutputPort()) + plane_mapper.ScalarVisibilityOn() + plane_mapper.SetScalarRange(data_scal.GetScalarRange()) + plane_mapper.SetLookupTable(clut) + plane_actor = vtk.vtkActor() + plane_actor.SetMapper(plane_mapper) + + #clip plane + plane_clip = vtk.vtkClipDataSet() + plane_clip.SetInputData(data_scal) + plane_clip.SetInputArrayToProcess( + 0, 0, 0, + vtk.vtkDataObject().FIELD_ASSOCIATION_POINTS, + vtk.vtkDataSetAttributes().SCALARS + ) + plane_clip.SetClipFunction(plane) + plane_clip.InsideOutOn() + clip_mapper = vtk.vtkDataSetMapper() + clip_mapper.SetInputConnection(plane_cut.GetOutputPort()) + clip_mapper.ScalarVisibilityOn() + clip_mapper.SetScalarRange(data_scal.GetScalarRange()) + clip_mapper.SetLookupTable(clut) + clip_actor = vtk.vtkActor() + clip_actor.SetMapper(clip_mapper) + + # Bounding box + outline_data = vtk.vtkOutlineFilter() + outline_data.SetInputData(data_test) + outline_mapper = vtk.vtkPolyDataMapper() + outline_mapper.SetInputConnection(outline_data.GetOutputPort()) + outline = vtk.vtkActor() + outline.SetMapper(outline_mapper) + outline.GetProperty().SetColor(0, 0, 0) + outline.GetProperty().SetLineWidth(5.0) + outline.SetPosition(-xmid, -ymid, -zmid) + + # Text + txt = vtk.vtkTextActor() + txt.SetInput('T = TEMP K') + txtprop = txt.GetTextProperty() + txtprop.SetFontFamilyToArial() + txtprop.SetFontSize(36) + txtprop.SetColor(0, 0, 0) + txt.SetDisplayPosition(30, 550) + + # Reposition the camera + ren.GetActiveCamera().Azimuth(0) + ren.GetActiveCamera().Elevation(0) + ren.GetActiveCamera().ParallelProjectionOn() + ren.GetActiveCamera().SetFocalPoint(0, 0, 0) + ren.GetActiveCamera().SetViewUp(0, 1, 0) + ren.GetActiveCamera().SetParallelScale(0.55 * ymax) + _length = max(xmax - xmin, zmax - zmin) / 2 + _height = _length / 0.26795 * 1.1 + + ren.GetActiveCamera().SetPosition(0, 0, _height) + ren.AddActor(vector) + + # Render scene + iren = vtk.vtkRenderWindowInteractor() + istyle = vtk.vtkInteractorStyleTrackballCamera() + iren.SetInteractorStyle(istyle) + iren.SetRenderWindow(ren_win) + iren.Initialize() + ren_win.Render() + screenshot(ren_win, temperature, field, plot_dir) + + +def plot_configuration(node: orm.Node, plot_dir: Union[str, os.PathLike]): + """Collect and plot the results""" + + for _sub_node in node.called_descendants: + + _index_temperature = node.inputs.temperatures.get_list().index(_sub_node.inputs.inpsd_dict.get_dict()['temp']) + _index_field = node.inputs.external_fields.get_list().index(_sub_node.inputs.inpsd_dict.get_dict()['hfield']) + + coord_array = _sub_node.get_node_by_label('coord') + + # Size of system for plot + max_number_atoms = 1000000 + mom_states = _sub_node.get_outgoing().get_node_by_label('mom_states_traj') + points, number_atoms = read_atoms(coord_array, max_number_atoms) + vectors, colors = read_vectors_data(mom_states, number_atoms) + if node.inputs.plot_individual.value > int(0): + plot_realspace(points, vectors, colors, _index_temperature, _index_field, plot_dir) + + +def plot_result(node: orm.Node, plot_dir: Union[str, os.PathLike], curie_temperature: Union[float, int, None] = None): + """ + Basic plot function + Note that we use normalized axis like T/Tc in all figures + """ + + fig_mag = go.Figure() + fig_ene = go.Figure() + + _norm_factor = curie_temperature if curie_temperature is not None else 1 + + for _sub_node in node.called_descendants: + data = _sub_node.outputs.temperature_output.get_dict() + fig_mag.add_trace( + go.Scatter( + np.array(data['temperature']) / _norm_factor, + np.array(data['magnetization']) / np.array(data['magnetization'][0]), + name=f'N={_sub_node.inputs.inpsd_dict.get_dict()["necll"].split()[0]}', + ) + ) + fig_ene.add_trace( + go.Scatter( + np.array(data['temperature']) / _norm_factor, + np.array(data['energy']), + name=f'N={_sub_node.inputs.inpsd_dict.get_dict()["necll"].split()[0]}', + ) + ) + + fig_mag.update_layout( + xaxis_title='T/Tc', + yaxis_title='M/Msat', + ) + fig_ene.update_layout( + xaxis_title='T/Tc', + yaxis_title='Energy (mRy)', + ) + + fig_mag.write_image(f'{plot_dir}/temperature-magnetization.png') + fig_ene.write_image(f'{plot_dir}/temperature-energy.png') + + +def plot_ams(node: orm.Node, plot_dir: Union[str, os.PathLike]): # pylint: disable=too-many-locals + """Process results and basic plot""" + + ams_plot_var_out = node.outputs.ams_plot_var.get_dict() + timestep = ams_plot_var_out['timestep'] + sc_step = ams_plot_var_out['sc_step'] + sqw_x = ams_plot_var_out['sqw_x'] + sqw_y = ams_plot_var_out['sqw_y'] + ams_t = ams_plot_var_out['ams'] + axidx_abs = ams_plot_var_out['axidx_abs'] + ams_dist_col = ams_plot_var_out['ams_dist_col'] + axlab = ams_plot_var_out['axlab'] + _ = plt.figure(figsize=[8, 5]) + axes = plt.subplot(111) + hbar = float(4.135667662e-15) + emax = 0.5 * float(hbar) / (float(timestep) * float(sc_step)) * 1e3 + #print('emax_sqw=',emax) + sqw_x = np.array(sqw_x) + sqw_y = np.array(sqw_y) + #ams = np.array(ams) + sqw_temp = (sqw_x**2.0 + sqw_y**2.0)**0.5 + sqw_temp[:, 0] = sqw_temp[:, 0] / 100.0 + sqw_temp = sqw_temp.T / sqw_temp.T.max(axis=0) + sqw_temp = ndimage.gaussian_filter1d( + sqw_temp, + sigma=1, + axis=1, + mode='constant', + ) + sqw_temp = ndimage.gaussian_filter1d( + sqw_temp, + sigma=5, + axis=0, + mode='reflect', + ) + plt.imshow( + sqw_temp, + cmap=plt.get_cmap('jet'), + interpolation='antialiased', + origin='lower', + extent=[axidx_abs[0], axidx_abs[-1], 0, emax], + ) # Note here, we could choose color bar here. + ams_t = np.array(ams_t) + plt.plot( + ams_t[:, 0] / ams_t[-1, 0] * axidx_abs[-1], + ams_t[:, 1:ams_dist_col], + 'white', + lw=1, + ) + plt.colorbar() + plt.xticks(axidx_abs, axlab) + plt.xlabel('q') + plt.ylabel('Energy (meV)') + plt.autoscale(tight=False) + axes.set_aspect('auto') + plt.grid(b=True, which='major', axis='x') + plt.savefig(f'{plot_dir}/AMS.png') + + +def plot_average_specific_heat(node: orm.Node, plot_dir: Union[str, os.PathLike]): + """Plot the average specific heat""" + + #plot the line for one B + + _energy = [] + _bfield = [] + _temperature = [] + for _sub_node in node.called_descendants: + _energy.append(_sub_node.outputs.cumulats.get_dict()['energy']) + _temperature.append(_sub_node.inputs.inpsd_dict.get_dict()['temp']) + _bfield.append(_sub_node.inputs.inpsd_dict.get_dict()['hfield']) + + _bfield_unique = np.unique(_bfield)[0] + + fig = go.Figure() + + heat_map = [] + + for _field in _bfield_unique: + _temperature_curr = np.asarray(_temperature)[np.where(np.asarray(_bfield) == _field)[0]] + _energy_curr = np.asarray(_energy[np.where(np.asarray(_bfield) == _field)[0]]) + + de_dt = np.gradient(np.array(_energy_curr)) / np.gradient(np.array(_temperature_curr)) + heat_map.append(de_dt.tolist()) + de_dt = (de_dt / de_dt[0]).tolist() + + fig.add_trace(go.Scatter(_temperature_curr, de_dt, name=f'B={_field}')) + fig.update_layout(xaxis_title='Temperature', yaxis_title='Cv') + fig.write_image(f'{plot_dir}/CV_T.png') + fig_heatmap = go.Figure(data=go.Heatmap( + z=heat_map, + x=_temperature, + y=_bfield, + zsmooth='best', + )) + + fig_heatmap.update_xaxes(range=[min(_temperature), max(_temperature)]) + fig_heatmap.update_yaxes(range=[int(min(_bfield)), int(max(_bfield))]) + fig_heatmap.update_traces(colorscale='Jet', selector=dict(type='heatmap')) + fig_heatmap.write_image(f'{plot_dir}/average_specific_heat.png') + + +def plot_average_magnetic_moment(node: orm.Node, plot_dir: Union[str, os.PathLike]): + """Plot the average magnetic moment""" + heat_map = [] + _temperature = [] + _bfield = [] + for _sub_node in node.called_descendats: + _temperature.append(_sub_node.inputs.inpsd_dict.get_dict()['temp']) + _bfield.append(_sub_node.inputs.inpsd_dcit.get_dict()['hfield']) + heat_map.append(_sub_node.outputs.cumulants.get_dict()['magnetization']) + + fig = go.Figure(data=go.Heatmap( + z=heat_map, + x=_temperature, + y=_bfield, + zsmooth='best', + )) + fig.update_xaxes(range=[min(_temperature), max(_temperature)]) + fig.update_yaxes(range=[int(min(_bfield)), int(max(_bfield))]) + fig.update_traces(colorscale='Jet', selector=dict(type='heatmap')) + fig.write_image(f'{plot_dir}/average_magnetic_moment.png') + + +def plot_pd( + plot_dir, + heat_map, + x_label_list, + y_label_list, + plot_name, + xlabel, + ylabel, +): + #pylint: disable=too-many-arguments + """Plotting the phase diagram""" + fig = go.Figure(data=go.Heatmap( + z=heat_map, + x=x_label_list, + y=y_label_list, + zsmooth='best', + )) + fig.update_xaxes(range=[int(float(x_label_list[0])), int(float(x_label_list[-1]))]) + fig.update_yaxes(range=[int(float(y_label_list[0])), int(float(y_label_list[-1]))]) + fig.update_traces(colorscale='Jet', selector=dict(type='heatmap')) + fig.update_layout(yaxis={'title': f'{xlabel}', 'tickangle': -90}, xaxis={'title': f'{ylabel}'}) + fig.write_image(f"{plot_dir}/{plot_name.replace(' ', '_')}.png") + + +def plot_line( + x_data, + y_data, + line_name_list, + x_label, + y_label, + plot_path, + plot_name, + leg_name_list, +): + #pylint: disable=too-many-arguments + """Line plot""" + #Since we need all lines in one plot here x and y should be a dict with the line name on that + plt.figure() + _, axes = plt.subplots() + for index, _entry in enumerate(line_name_list): + axes.plot( + x_data, + y_data[_entry], + label=f'{leg_name_list[index]}', + ) + axes.legend() + axes.set_xlabel(x_label) + axes.set_ylabel(y_label) + plt.savefig(f'{plot_path}/{plot_name}.png') + plt.close() + + +def plot_mag_phase_diagram(node: orm.Node): + """Plot the magnetization phase diagram""" + y_label_list = [] + for i in node.inputs.external_fields.get_list(): + y_label_list.append(float(max(np.array(i.split())))) + + for _, cell_size in enumerate(node.inputs.cell_size): + pd_dict = node.outputs.thermal_dynamic_output.get_dict()[f'{cell_size}'] + pd_for_plot = [] + for i in pd_dict.keys(): + pd_for_plot.append(pd_dict[i]['magnetization']) + plot_pd( + node.inputs.plot_dir.value, pd_for_plot, node.inputs.temperatures.get_list(), y_label_list, + ('M_T' + str(cell_size)), 'B', 'T' + ) + + +def plot_mag_temp(node: orm.Node): + """Plot the magnetization as a function of temperature""" + line_name_list = [] + leg_name_list = [] + line_for_plot = {} + for _, cell_size in enumerate(node.inputs.cell_size): + line_dict = node.outputs.thermal_dynamic_output.get_dict()[f'{cell_size}'] + leg_name_list.append(('Cell:' + cell_size.replace(' ', '_'))) + index = 0 + for i in line_dict.keys(): + line_for_plot[f"{cell_size.replace(' ', '_')}_{index}"] = line_dict[i]['magnetization'] + line_name_list.append(f"{cell_size.replace(' ', '_')}_{index}") + index = index + 1 + plot_line( + node.inputs.temperatures.get_list(), + line_for_plot, + line_name_list, + 'T', + 'M', + node.inputs.plot_dir.value, + 'M_T', + leg_name_list, + ) + + +def plot_cev_phase_diagram(node: orm.Node): + """Check if the specific heat phase diagram should be produced""" + y_label_list = [] + for i in node.inputs.external_fields.get_list(): + y_label_list.append(float(max(np.array(i.split())))) + + for _, cell_size in enumerate(node.inputs.cell_size): + pd_dict = node.outputs.thermal_dynamic_output.get_dict()[f'{cell_size}'] + pd_for_plot = [] + for i in pd_dict.keys(): + pd_for_plot.append(pd_dict[i]['specific_heat']) + plot_pd( + node.inputs.plot_dir.value, pd_for_plot, node.inputs.temperatures.get_list(), y_label_list, + ('specific_heat_T' + str(cell_size)), 'B', 'T' + ) + + +def plot_specific_heat_temp(node: orm.Node): + """Plot the specific head as as function of temperature""" + line_name_list = [] + leg_name_list = [] + line_for_plot = {} + for _, cell_size in enumerate(node.inputs.cell_size): + line_dict = node.outputs.thermal_dynamic_output.get_dict()[f'{cell_size}'] + leg_name_list.append(('Cell:' + cell_size.replace(' ', '_'))) + index = 0 + for i in line_dict.keys(): + line_for_plot[f"{cell_size.replace(' ', '_')}_{index}"] = line_dict[i]['specific_heat'] + line_name_list.append(f"{cell_size.replace(' ', '_')}_{index}") + index = index + 1 + plot_line( + node.inputs.temperatures.get_list(), line_for_plot, line_name_list, 'T', 'specific_heat', + node.inputs.plot_dir.value, 'specific_heat_T', leg_name_list + ) + + +def plot_sus_phase_diagram(node: orm.Node): + """Plot the magnetic susceptibility phase diagram""" + y_label_list = [] + for i in node.inputs.external_fields.get_list(): + y_label_list.append(float(max(np.array(i.split())))) + + for _, cell_size in enumerate(node.inputs.cell_size): + pd_dict = node.outputs.thermal_dynamic_output.get_dict()[f'{cell_size}'] + pd_for_plot = [] + for i in pd_dict.keys(): + pd_for_plot.append(pd_dict[i]['susceptibility']) + plot_pd( + node.inputs.plot_dir.value, pd_for_plot, node.inputs.temperatures.get_list(), y_label_list, + ('susceptibility_T' + str(cell_size)), 'B', 'T' + ) + + +def plot_susceptibility_temp(node: orm.Node): + """Plot the susceptibility as a function of temperature""" + line_name_list = [] + leg_name_list = [] + line_for_plot = {} + for _, cell_size in enumerate(node.inputs.cell_size): + line_dict = node.outputs.thermal_dynamic_output.get_dict()[f'{cell_size}'] + leg_name_list.append(('Cell:' + cell_size.replace(' ', '_'))) + index = 0 + for i in line_dict.keys(): + line_for_plot[f"{cell_size.replace(' ', '_')}_{index}"] = line_dict[i]['susceptibility'] + line_name_list.append(f"{cell_size.replace(' ', '_')}_{index}") + index = index + 1 + plot_line( + node.inputs.temperatures.get_list(), line_for_plot, line_name_list, 'T', 'Susceptibility', + node.inputs.plot_dir.value, 'Susceptibility_T', leg_name_list + ) + + +def plot_free_e_phase_diagram(node: orm.Node): + """Plot the free energy phase diagram""" + y_label_list = [] + for i in node.inputs.external_fields.get_list(): + y_label_list.append(float(max(np.array(i.split())))) + + for _, cell_size in enumerate(node.inputs.cell_size): + pd_dict = node.outputs.thermal_dynamic_output.get_dict()[f'{cell_size}'] + pd_for_plot = [] + for i in pd_dict.keys(): + pd_for_plot.append(pd_dict[i]['free_e']) + plot_pd( + node.inputs.plot_dir.value, + pd_for_plot, + node.inputs.temperatures.get_list(), + y_label_list, + ('free_e_T' + str(cell_size)), + 'B', + 'T', + ) + + +def plot_free_e_temp(node: orm.Node): + """Plot the free energy as a function of temperature""" + line_name_list = [] + leg_name_list = [] + line_for_plot = {} + for _, cell_size in enumerate(node.inputs.cell_size): + line_dict = node.outputs.thermal_dynamic_output.get_dict()[f'{cell_size}'] + leg_name_list.append(('Cell:' + cell_size.replace(' ', '_'))) + index = 0 + for i in line_dict.keys(): + line_for_plot[f"{cell_size.replace(' ', '_')}_{index}"] = line_dict[i]['free_e'] + line_name_list.append(f"{cell_size.replace(' ', '_')}_{index}") + index = index + 1 + plot_line( + node.inputs.temperatures.get_list(), line_for_plot, line_name_list, 'T', 'free_e', node.inputs.plot_dir.value, + 'free_e_T', leg_name_list + ) + + +def plot_entropy_phase_diagram(node: orm.Node): + """Plot the entropy phase diagram""" + y_label_list = [] + for i in node.inputs.external_fields.get_list(): + y_label_list.append(float(max(np.array(i.split())))) + + for _, cell_size in enumerate(node.inputs.cell_size): + pd_dict = node.outputs.thermal_dynamic_output.get_dict()[f'{cell_size}'] + pd_for_plot = [] + for i in pd_dict.keys(): + pd_for_plot.append(pd_dict[i]['entropy']) + plot_pd( + node.inputs.plot_dir.value, pd_for_plot, node.inputs.temperatures.get_list(), y_label_list, + ('entropy_T' + str(cell_size)), 'B', 'T' + ) + + +def plot_entropy_temp(node: orm.Node): + """Plot the entropy vs temperature""" + line_name_list = [] + leg_name_list = [] + line_for_plot = {} + for _, cell_size in enumerate(node.inputs.cell_size): + line_dict = node.outputs.thermal_dynamic_output.get_dict()[f'{cell_size}'] + leg_name_list.append(('Cell:' + cell_size.replace(' ', '_'))) + index = 0 + for i in line_dict.keys(): + line_for_plot[f"{cell_size.replace(' ', '_')}_{index}"] = line_dict[i]['entropy'] + line_name_list.append(f"{cell_size.replace(' ', '_')}_{index}") + index = index + 1 + plot_line( + node.inputs.temperatures.get_list(), line_for_plot, line_name_list, 'T', 'Entropy', node.inputs.plot_dir.value, + 'Entropy_T', leg_name_list + ) + + +def plot_energy_phase_diagram(node: orm.Node): + """Plot the energy phase diagram""" + y_label_list = [] + for i in node.inputs.external_fields.get_list(): + y_label_list.append(float(max(np.array(i.split())))) + + for _, cell_size in enumerate(node.inputs.cell_size): + pd_dict = node.outputs.thermal_dynamic_output.get_dict()[f'{cell_size}'] + pd_for_plot = [] + for i in pd_dict.keys(): + pd_for_plot.append(pd_dict[i]['energy']) + plot_pd( + node.inputs.plot_dir.value, pd_for_plot, node.inputs.temperatures.get_list(), y_label_list, + ('energy_T' + str(cell_size)), 'B', 'T' + ) + + +def plot_energy_temp(node: orm.Node): + """Plot the energy vs temperature""" + line_name_list = [] + leg_name_list = [] + line_for_plot = {} + for _, cell_size in enumerate(node.inputs.cell_size): + line_dict = node.outputs.thermal_dynamic_output.get_dict()[f'{cell_size}'] + leg_name_list.append(('Cell:' + cell_size.replace(' ', '_'))) + index = 0 + for i in line_dict.keys(): + line_for_plot[f"{cell_size.replace(' ', '_')}_{index}"] = line_dict[i]['energy'] + line_name_list.append(f"{cell_size.replace(' ', '_')}_{index}") + index = index + 1 + plot_line( + node.inputs.temperatures.get_list(), line_for_plot, line_name_list, 'T', 'energy', node.inputs.plot_dir.value, + 'energy_T', leg_name_list + ) + + +def plot_dudt_phase_diagram(node: orm.Node): + """Plot the variation of the energy phase diagram""" + y_label_list = [] + for i in node.inputs.external_fields.get_list(): + y_label_list.append(float(max(np.array(i.split())))) + + for _, cell_size in enumerate(node.inputs.cell_size): + pd_dict = node.outputs.thermal_dynamic_output.get_dict()[f'{cell_size}'] + pd_for_plot = [] + for i in pd_dict.keys(): + pd_for_plot.append(pd_dict[i]['dudt']) + plot_pd( + node.inputs.plot_dir.value, pd_for_plot, node.inputs.temperatures.get_list(), y_label_list, + ('dudt_T' + str(cell_size)), 'B', 'T' + ) + + +def plot_dudt_temp(node: orm.Node): + """Plot the variation of the energy vs temperature""" + line_name_list = [] + leg_name_list = [] + line_for_plot = {} + for _, cell_size in enumerate(node.inputs.cell_size): + line_dict = node.outputs.thermal_dynamic_output.get_dict()[f'{cell_size}'] + leg_name_list.append(('Cell:' + cell_size.replace(' ', '_'))) + index = 0 + for i in line_dict.keys(): + line_for_plot[f"{cell_size.replace(' ', '_')}_{index}"] = line_dict[i]['dudt'] + line_name_list.append(f"{cell_size.replace(' ', '_')}_{index}") + index = index + 1 + plot_line( + node.inputs.temperatures.get_list(), + line_for_plot, + line_name_list, + 'T', + 'dudt', + node.inputs.plot_dir.value, + 'dudt_T', + leg_name_list, + ) + + +def plot_binder_cumu_phase_diagram(node: orm.Node): + """Plot the binder cumulant phase diagram""" + y_label_list = [] + for i in node.inputs.external_fields.get_list(): + y_label_list.append(float(max(np.array(i.split())))) + + for _, cell_size in enumerate(node.inputs.cell_size): + pd_dict = node.outputs.thermal_dynamic_output.get_dict()[f'{cell_size}'] + pd_for_plot = [] + for i in pd_dict.keys(): + pd_for_plot.append(pd_dict[i]['binder_cumulant']) + plot_pd( + node.inputs.plot_dir.value, pd_for_plot, node.inputs.temperatures.get_list(), y_label_list, + ('binder_cumulant_T' + str(cell_size)), 'B', 'T' + ) + + +def plot_binder_cumulant_temp(node: orm.Node): + """Plot the binder cumulant as as function of temperature""" + line_name_list = [] + leg_name_list = [] + line_for_plot = {} + for _, cell_size in enumerate(node.inputs.cell_size): + line_dict = node.outputs.thermal_dynamic_output.get_dict()[f'{cell_size}'] + leg_name_list.append(('Cell:' + cell_size.replace(' ', '_'))) + index = 0 + for i in line_dict.keys(): + line_for_plot[f"{cell_size.replace(' ', '_')}_{index}"] = line_dict[i]['binder_cumulant'] + line_name_list.append(f"{cell_size.replace(' ', '_')}_{index}") + index = index + 1 + plot_line( + node.inputs.temperatures.get_list(), + line_for_plot, + line_name_list, + 'T', + 'binder_cumulant', + node.inputs.plot_dir.value, + 'binder_cumulant_T', + leg_name_list, + ) diff --git a/aiida_uppasd/workflows/data_generating.py b/aiida_uppasd/workflows/data_generating.py index bf0de1b..24721d8 100644 --- a/aiida_uppasd/workflows/data_generating.py +++ b/aiida_uppasd/workflows/data_generating.py @@ -7,16 +7,12 @@ Workchain demo generating data for Machine Learning """ -from PIL import Image -import matplotlib.pyplot as plt import numpy as np -import plotly.graph_objects as go from aiida.engine import ToContext, WorkChain, calcfunction, if_ from aiida.orm import ArrayData, Dict, Int, List, Str from aiida_uppasd.workflows.base_restart import ASDBaseRestartWorkChain -from aiida_uppasd.workflows.skyrmions_pd_and_graph import plot_realspace, read_atoms, read_vectors_data @calcfunction @@ -108,18 +104,7 @@ def define(cls, spec): ) spec.outline( cls.submits, - if_(cls.check_for_plot_sk)( - cls.results_and_plot, - cls.combined_for_pd, - ), - if_(cls.check_for_plot_sk_number)( - cls.plot_skynumber, - cls.plot_skynumber_avg, - cls.plot_skynumber_number_heat_map, - cls.plot_sknumber_ave_heatmap, - ), - if_(cls.check_plot_ave_magnetic_mom)(cls.plot_average_magnetic_moment), - if_(cls.check_plot_ave_specific_heat)(cls.plot_average_specific_heat), + if_(cls.check_for_plot_sk)(cls.results_and_plot,), ) def check_for_plot_sk(self): @@ -131,33 +116,6 @@ def check_for_plot_sk(self): except BaseException: return False - def check_for_plot_sk_number(self): - """Check if the skyrmion number should be plotted""" - try: - if self.inputs.sk_number_plot.value > int(0): - return True - return False - except BaseException: - return False - - def check_plot_ave_magnetic_mom(self): - """Check fi the average magnetization should be plotted""" - try: - if self.inputs.average_magnetic_moment_plot.value > int(0): - return True - return False - except BaseException: - return False - - def check_plot_ave_specific_heat(self): - """Check if the average specific heat should be plotted""" - try: - if self.inputs.average_specific_heat_plot.value > int(0): - return True - return False - except BaseException: - return False - def generate_inputs(self): """Generate the inputs for the calculation""" inputs = self.exposed_inputs( @@ -197,23 +155,6 @@ def results_and_plot(self): sm_list = [] sam_list = [] for index_field, _ in enumerate(self.inputs.external_fields): - coord_array = self.ctx[f'T{index_temperature}_B{index_field}'].get_outgoing().get_node_by_label('coord') - - # Size of system for plot - max_number_atoms = 1000000 - mom_states = self.ctx[f'T{index_temperature}_B{index_field}'].get_outgoing( - ).get_node_by_label('mom_states_traj') - points, number_atoms = read_atoms(coord_array, max_number_atoms) - vectors, colors = read_vectors_data(mom_states, number_atoms) - plot_dir = self.inputs.plot_dir.value - plot_realspace( - points, - vectors, - colors, - index_temperature, - index_field, - plot_dir, - ) sk_data = self.ctx[f'T{index_temperature}_B{index_field}'].get_outgoing( ).get_node_by_label('sk_num_out') @@ -226,204 +167,3 @@ def results_and_plot(self): sk_num_for_plot = store_sknum_data(List(list=sky_num_out), List(list=sky_avg_num_out)) self.out('sk_num_for_plot', sk_num_for_plot) - - def combined_for_pd(self): - """Combine the data for the phase diagram""" - skyrmions_singleplot = [] - plot_dir = self.inputs.plot_dir.value - for index_temperature, _ in enumerate(self.inputs.temperatures): - for index_field, _ in enumerate(self.inputs.external_fields): - skyrmions_singleplot.append( - Image.open(f'{plot_dir}/T{index_temperature}B{index_field}.png').crop((200, 0, 848, 640)) - ) #we need to crop the white edges - phase_diagram = Image.new('RGB', (648 * len(self.inputs.temperatures), 640 * len(self.inputs.external_fields))) - for index_temperature, _ in enumerate(self.inputs.temperatures): - for index_field, _ in enumerate(self.inputs.external_fields): - phase_diagram.paste( - skyrmions_singleplot[len(self.inputs.external_fields) * index_temperature + index_field], - (0 + 648 * index_temperature, 0 + 640 * index_field) - ) - phase_diagram.save(f'{plot_dir}/PhaseDiagram.png', quality=100) - - def plot_skynumber(self): - """Plot the skyrmion number""" - plot_dir = self.inputs.plot_dir.value - _, axes = plt.subplots() - for index_field, external_field in enumerate(self.inputs.external_fields): - sk_number_list = [] - for index_temperature, _ in enumerate(self.inputs.temperatures): - sk_number_list.append( - self.ctx[f'T{index_temperature}_B{index_field}'].get_outgoing().get_node_by_label('sk_num_out'). - get_array('sk_number_data')[0] - ) - axes.plot(self.inputs.temperatures.get_list(), sk_number_list, label=f'B: {external_field[-6:]} T') - axes.legend(fontsize='xx-small') - axes.set_ylabel('Skyrmion number') - axes.set_xlabel('Temperature') - plt.savefig(f'{plot_dir}/sk_number.png') - - def plot_skynumber_avg(self): - """Plot the average skyrmion number""" - plot_dir = self.inputs.plot_dir.value - _, axes = plt.subplots() - for index_field, external_field in enumerate(self.inputs.external_fields): - sk_number_list = [] - for index_temperature, _ in enumerate(self.inputs.temperatures): - sk_number_list.append( - self.ctx[f'T{index_temperature}_B{index_field}'].get_outgoing().get_node_by_label('sk_num_out'). - get_array('sk_number_data')[1] - ) - axes.plot( - self.inputs.temperatures.get_list(), sk_number_list, label=f'B: {external_field[-6:]} T' - ) #this should be changed with different B - axes.legend(fontsize='xx-small') - axes.set_ylabel('Skyrmion number') - axes.set_xlabel('Temperature') - plt.savefig(f'{plot_dir}/sk_number_avg.png') - - def plot_skynumber_number_heat_map(self): - """Plot the average skyrmion number heat map""" - plot_dir = self.inputs.plot_dir.value - heat_map = [] - for index_field, _ in enumerate(self.inputs.external_fields): - sk_number_list = [] - for index_temperature, _ in enumerate(self.inputs.temperatures): - sk_number_list.append( - self.ctx[f'T{index_temperature}_B{index_field}'].get_outgoing().get_node_by_label('sk_num_out'). - get_array('sk_number_data')[0] - ) - heat_map.append(sk_number_list) - y_orginal = self.inputs.external_fields.get_list() - y_for_plot = [] - for i in y_orginal: - y_for_plot.append(float(i.split()[-1])) - - fig = go.Figure( - data=go.Heatmap( - z=heat_map, - x=self.inputs.temperatures.get_list(), - y=y_for_plot, - zsmooth='best', - ) - ) - - fig.update_xaxes(range=[self.inputs.temperatures.get_list()[0], self.inputs.temperatures.get_list()[-1]]) - fig.update_yaxes(range=[y_for_plot[0], y_for_plot[-1]]) - fig.update_traces(colorscale='Jet', selector=dict(type='heatmap')) - fig.write_image(f'{plot_dir}/sk_number_heatmap.png') - - def plot_sknumber_ave_heatmap(self): - """Plot the skyrmion number average heat map""" - plot_dir = self.inputs.plot_dir.value - heat_map = [] - for index_field, _ in enumerate(self.inputs.external_fields): - sk_number_list = [] - for index_temperature, _ in enumerate(self.inputs.temperatures): - sk_number_list.append( - self.ctx[f'T{index_temperature}_B{index_field}'].get_outgoing().get_node_by_label('sk_num_out'). - get_array('sk_number_data')[1] - ) - heat_map.append(sk_number_list) - y_orginal = self.inputs.external_fields.get_list() - y_for_plot = [] - for i in y_orginal: - y_for_plot.append(float(i.split()[-1])) - - fig = go.Figure( - data=go.Heatmap( - z=heat_map, - x=self.inputs.temperatures.get_list(), - y=y_for_plot, - zsmooth='best', - ) - ) - - fig.update_xaxes(range=[self.inputs.temperatures.get_list()[0], self.inputs.temperatures.get_list()[-1]]) - fig.update_yaxes(range=[y_for_plot[0], y_for_plot[-1]]) - fig.update_traces(colorscale='Jet', selector=dict(type='heatmap')) - fig.write_image(f'{plot_dir}/sk_number_heatmap_avg.png') - - def plot_average_magnetic_moment(self): - """Plot the average magnetic moment""" - plot_dir = self.inputs.plot_dir.value - heat_map = [] - for index_field, _ in enumerate(self.inputs.external_fields): - average_magnetic_moment = [] - for index_temperature, _ in enumerate(self.inputs.temperatures): - average_magnetic_moment.append( - self.ctx[f'T{index_temperature}_B{index_field}'].get_outgoing().get_node_by_label('cumulants') - ['magnetization'] - ) - heat_map.append(average_magnetic_moment) - y_orginal = self.inputs.external_fields.get_list() - y_for_plot = [] - for i in y_orginal: - y_for_plot.append(float(i.split()[-1])) - - fig = go.Figure( - data=go.Heatmap( - z=heat_map, - x=self.inputs.temperatures.get_list(), - y=y_for_plot, - zsmooth='best', - ) - ) - fig.update_xaxes(range=[self.inputs.temperatures.get_list()[0], self.inputs.temperatures.get_list()[-1]]) - fig.update_yaxes(range=[int(y_for_plot[0]), int(y_for_plot[-1])]) - fig.update_traces(colorscale='Jet', selector=dict(type='heatmap')) - fig.write_image(f'{plot_dir}/average_magnetic_moment.png') - - def plot_average_specific_heat(self): - """Plot the average specific heat""" - plot_dir = self.inputs.plot_dir.value - - if len(self.inputs.external_fields.get_list()) == 1: - #plot the line for one B - for index_field, _ in enumerate(self.inputs.external_fields): - energy_list = [] - temperature_list = self.inputs.temperatures.get_list() - for index_temperature, _ in enumerate(self.inputs.temperatures): - energy_list.append( - self.ctx[f'T{index_temperature}_B{index_field}'].get_outgoing().get_node_by_label('cumulants') - ['energy'] - ) - de_dt = (np.gradient(np.array(energy_list)) / np.gradient(np.array(temperature_list))) - de_dt = (de_dt / de_dt[0]).tolist() - fig, axes = plt.subplots() - axes.plot(self.inputs.temperatures.get_list(), de_dt) - axes.legend(fontsize='xx-small') - axes.set_ylabel('C') - axes.set_xlabel('Temperature') - plt.savefig(f'{plot_dir}/CV_T.png') - - else: - heat_map = [] - for index_field, _ in enumerate(self.inputs.external_fields): - energy_list = [] - temperature_list = self.inputs.temperatures.get_list() - - for index_temperature, _ in enumerate(self.inputs.temperatures): - energy_list.append( - self.ctx[f'T{index_temperature}_B{index_field}'].get_outgoing().get_node_by_label('cumulants') - ['energy'] - ) - de_dt = (np.gradient(np.array(energy_list)) / np.gradient(np.array(temperature_list))).tolist() - heat_map.append(de_dt) - y_orginal = self.inputs.external_fields.get_list() - y_for_plot = [] - for i in y_orginal: - y_for_plot.append(float(i.split()[-1])) - - fig = go.Figure( - data=go.Heatmap( - z=heat_map, - x=self.inputs.temperatures.get_list(), - y=y_for_plot, - zsmooth='best', - ) - ) - - fig.update_xaxes(range=[self.inputs.temperatures.get_list()[0], self.inputs.temperatures.get_list()[-1]]) - fig.update_yaxes(range=[int(y_for_plot[0]), int(y_for_plot[-1])]) - fig.update_traces(colorscale='Jet', selector=dict(type='heatmap')) - fig.write_image(f'{plot_dir}/average_specific_heat.png') diff --git a/aiida_uppasd/workflows/magnon_spectra.py b/aiida_uppasd/workflows/magnon_spectra.py index 8ffd729..f592639 100644 --- a/aiida_uppasd/workflows/magnon_spectra.py +++ b/aiida_uppasd/workflows/magnon_spectra.py @@ -17,11 +17,6 @@ which @Anders bergman """ -import matplotlib.pyplot as plt -# -*- coding: utf-8 -*- -import numpy as np -from scipy import ndimage - from aiida import orm from aiida.engine import ToContext, WorkChain @@ -71,6 +66,12 @@ def define(cls, spec): required=True, help='Adiabatic magnon spectrum', ) + spec.output( + 'ams_plot_var', + valid_type=orm.Dict, + required=True, + help='Adiabatic magnon spectrum plotting variables', + ) spec.outline( cls.submits, @@ -118,62 +119,9 @@ def submits(self): calculations['AMS'] = future return ToContext(**calculations) - def results_and_plot(self): # pylint: disable=too-many-locals + def results(self): # pylint: disable=too-many-locals """Process results and basic plot""" ams_plot_var_out = self.ctx['AMS'].get_outgoing().get_node_by_label('AMS_plot_var') - timestep = ams_plot_var_out.get_dict()['timestep'] - sc_step = ams_plot_var_out.get_dict()['sc_step'] - sqw_x = ams_plot_var_out.get_dict()['sqw_x'] - sqw_y = ams_plot_var_out.get_dict()['sqw_y'] - ams_t = ams_plot_var_out.get_dict()['ams'] - axidx_abs = ams_plot_var_out.get_dict()['axidx_abs'] - ams_dist_col = ams_plot_var_out.get_dict()['ams_dist_col'] - axlab = ams_plot_var_out.get_dict()['axlab'] - plot_dir = self.inputs.plot_dir.value - _ = plt.figure(figsize=[8, 5]) - axes = plt.subplot(111) - hbar = float(4.135667662e-15) - emax = 0.5 * float(hbar) / (float(timestep) * float(sc_step)) * 1e3 - #print('emax_sqw=',emax) - sqw_x = np.array(sqw_x) - sqw_y = np.array(sqw_y) - #ams = np.array(ams) - sqw_temp = (sqw_x**2.0 + sqw_y**2.0)**0.5 - sqw_temp[:, 0] = sqw_temp[:, 0] / 100.0 - sqw_temp = sqw_temp.T / sqw_temp.T.max(axis=0) - sqw_temp = ndimage.gaussian_filter1d( - sqw_temp, - sigma=1, - axis=1, - mode='constant', - ) - sqw_temp = ndimage.gaussian_filter1d( - sqw_temp, - sigma=5, - axis=0, - mode='reflect', - ) - plt.imshow( - sqw_temp, - cmap=plt.get_cmap('jet'), - interpolation='antialiased', - origin='lower', - extent=[axidx_abs[0], axidx_abs[-1], 0, emax], - ) # Note here, we could choose color bar here. - ams_t = np.array(ams_t) - plt.plot( - ams_t[:, 0] / ams_t[-1, 0] * axidx_abs[-1], - ams_t[:, 1:ams_dist_col], - 'white', - lw=1, - ) - plt.colorbar() - plt.xticks(axidx_abs, axlab) - plt.xlabel('q') - plt.ylabel('Energy (meV)') - plt.autoscale(tight=False) - axes.set_aspect('auto') - plt.grid(b=True, which='major', axis='x') - plt.savefig(f'{plot_dir}/AMS.png') ams = self.ctx['AMS'].get_outgoing().get_node_by_label('ams') self.out('ams', ams) + self.out('ams_plot_var', ams_plot_var_out) diff --git a/aiida_uppasd/workflows/skyrmions_pd_and_graph.py b/aiida_uppasd/workflows/skyrmions_pd_and_graph.py index 4e1d777..2ffc59c 100644 --- a/aiida_uppasd/workflows/skyrmions_pd_and_graph.py +++ b/aiida_uppasd/workflows/skyrmions_pd_and_graph.py @@ -1,4 +1,3 @@ -#!/usr/bin/env python3 # -*- coding: utf-8 -*- """ Created on Tue Feb 22 13:03:40 2022 @@ -17,11 +16,7 @@ which @Anders bergman """ -from PIL import Image -import matplotlib.pyplot as plt import numpy as np -import plotly.graph_objects as go -import vtk from aiida.engine import ToContext, WorkChain, calcfunction, if_ from aiida.orm import ArrayData, Dict, Int, List, Str @@ -29,294 +24,6 @@ from aiida_uppasd.workflows.base_restart import ASDBaseRestartWorkChain -# Read Location of Atoms -def read_atoms(coord_array_data, max_number_atoms): - """Transform the atomic positions to VTK structure""" - - points = vtk.vtkPoints() - number_atoms = 0 - # Read ahead - - # Read all data - for data in coord_array_data.get_array('coord'): - if number_atoms <= max_number_atoms: - coord_x, coord_y, coord_z = float(data[1]), float(data[2]), float(data[3]) - points.InsertPoint(number_atoms, coord_x, coord_y, coord_z) - number_atoms = number_atoms + 1 - return points, number_atoms - - -# Read vectors -# We must know the time step and the number of atoms per time -def read_vectors_data(mom_states, number_atoms): - """Transform the magnetic moments to VTK structure""" - final_mom_states_x = mom_states.get_array('mom_states_x')[-number_atoms:] - final_mom_states_y = mom_states.get_array('mom_states_y')[-number_atoms:] - final_mom_states_z = mom_states.get_array('mom_states_z')[-number_atoms:] - # Create a Double array which represents the vectors - vectors = vtk.vtkFloatArray() - colors = vtk.vtkFloatArray() - - # Define number of elemnts - vectors.SetNumberOfComponents(3) - colors.SetNumberOfComponents(1) - - for index in range(number_atoms): - vec_x, vec_y, vec_z = float(final_mom_states_x[index]), float(final_mom_states_y[index] - ), float(final_mom_states_z[index]) - _color = (vec_z + 1.0) / 2.0 - vectors.InsertTuple3(index, vec_x, vec_y, vec_z) - colors.InsertValue(index, _color) - return vectors, colors - - -#---------------------------- -# screenshot code begins here -#---------------------------- -#A function that takes a renderwindow and saves its contents to a .png file -def screenshot(ren_win, temperature, field, plot_dir): - """Take a screenshot of the magnetic configuration and save it as a png""" - win2im = vtk.vtkWindowToImageFilter() - win2im.ReadFrontBufferOff() - win2im.SetInput(ren_win) - # - povexp = vtk.vtkPOVExporter() - povexp.SetRenderWindow(ren_win) - #povexp.SetInput(renWin) - ren_win.Render() - # ren_win.SetFileName('realspace_plot_T_{}_B_{}.pov'.format(T,B)) - # povexp.Write() - # - to_png = vtk.vtkPNGWriter() - to_png.SetFileName(f'{plot_dir}/T{temperature}B{field}.png') - to_png.SetInputConnection(win2im.GetOutputPort()) - to_png.Write() - return 0 - - -def plot_realspace(points, vectors, colors, temperature, field, plot_dir): #pylint: disable=too-many-statements,too-many-locals - """Plot the real space magnetic configuration""" - ren_win = vtk.vtkRenderWindow() - ren = vtk.vtkOpenGLRenderer() - ren_win.AddRenderer(ren) - - # Set color of backgroung - ren.SetBackground(1.0, 1.0, 1.0) - ren_win.SetSize(2096, 1280) - - data_test = vtk.vtkPolyData() - data_scal = vtk.vtkUnstructuredGrid() - - # Read atom positions - atom_data = points - data_test.SetPoints(atom_data) - data_scal.SetPoints(atom_data) - - # Read data for vectors - vecz, colz = (vectors, colors) - data_test.GetPointData().SetVectors(vecz) - data_test.GetPointData().SetScalars(colz) - data_scal.GetPointData().SetScalars(colz) - - # Create colortable for the coloring of the vectors - lut = vtk.vtkLookupTable() - for i in range(0, 128, 1): - lut.SetTableValue(i, (127.0 - i) / 127.0, i / 127.0, 0, 1) - for i in range(128, 256, 1): - lut.SetTableValue(i, 0, (256.0 - i) / 128.0, (i - 128.0) / 128, 1) - lut.SetTableRange(-1.0, 1.0) - lut.Build() - - # Set up atoms - ball = vtk.vtkSphereSource() - ball.SetRadius(1.00) - ball.SetThetaResolution(16) - ball.SetPhiResolution(16) - - balls = vtk.vtkGlyph3DMapper() - balls.SetInputData(data_test) - balls.SetSourceConnection(ball.GetOutputPort()) - balls.SetScaleFactor(0.15) - balls.SetScaleModeToNodata_scaling() - balls.SetLookupTable(lut) - balls.Update() - - atom = vtk.vtkLODActor() - atom.SetMapper(balls) - atom.GetProperty().SetOpacity(1.5) - xmin, xmax = atom.GetXRange() - ymin, ymax = atom.GetYRange() - zmin, zmax = atom.GetZRange() - xmid = (xmin + xmax) / 2 - ymid = (ymin + ymax) / 2 - zmid = (zmin + zmax) / 2 - atom.SetPosition(-xmid, -ymid, -zmid) - atom.GetProperty().SetSpecular(0.3) - atom.GetProperty().SetSpecularPower(60) - atom.GetProperty().SetAmbient(0.2) - atom.GetProperty().SetDiffuse(0.8) - - # Create vectors - arrow = vtk.vtkArrowSource() - arrow.SetTipRadius(0.25) - arrow.SetShaftRadius(0.15) - arrow.SetTipResolution(72) - arrow.SetShaftResolution(72) - - glyph = vtk.vtkGlyph3DMapper() - glyph.SetSourceConnection(arrow.GetOutputPort()) - glyph.SetInputData(data_test) - glyph.SetScaleFactor(2.00) - # Color the vectors according to magnitude - glyph.SetScaleModeToNodata_scaling() - glyph.SetLookupTable(lut) - glyph.SetColorModeToMapScalars() - glyph.Update() - - vector = vtk.vtkLODActor() - vector.SetMapper(glyph) - vector.SetPosition(-xmid, -ymid, -zmid) - - vector.GetProperty().SetSpecular(0.3) - vector.GetProperty().SetSpecularPower(60) - vector.GetProperty().SetAmbient(0.2) - vector.GetProperty().SetDiffuse(0.8) - vector.GetProperty().SetOpacity(1.0) - - # Scalar bar - scalar_bar = vtk.vtkScalarBarActor() - scalar_bar.SetLookupTable(lut) - scalar_bar.SetOrientationToHorizontal() - scalar_bar.SetNumberOfLabels(0) - scalar_bar.SetPosition(0.1, 0.05) - scalar_bar.SetWidth(0.85) - scalar_bar.SetHeight(0.3) - scalar_bar.GetLabelTextProperty().SetFontSize(8) - - #Depth sorted field - data_test = vtk.vtkDepthSortPolyData() - data_test.SetInputData(data_test) - data_test.SetCamera(ren.GetActiveCamera()) - data_test.SortScalarsOn() - data_test.Update() - - #cubes - cube = vtk.vtkCubeSource() - cubes = vtk.vtkGlyph3DMapper() - cubes.SetInputConnection(data_test.GetOutputPort()) - cubes.SetSourceConnection(cube.GetOutputPort()) - cubes.SetScaleModeToNodata_scaling() - cubes.SetScaleFactor(0.995) - cubes.SetLookupTable(lut) - cubes.SetColorModeToMapScalars() - cubes.ScalarVisibilityOn() - cubes.OrientOff() - cubes.Update() - cubes.SetLookupTable(lut) - - cube_actor = vtk.vtkActor() - cube_actor.SetMapper(cubes) - cube_actor.GetProperty().SetOpacity(0.05) - cube_actor.SetPosition(-xmid, -ymid, -zmid) - - #hedgehog - hhog = vtk.vtkHedgeHog() - hhog.SetInputData(data_test) - hhog.SetScaleFactor(5.0) - hhog_mapper = vtk.vtkPolyDataMapper() - hhog_mapper.SetInputConnection(hhog.GetOutputPort()) - hhog_mapper.SetLookupTable(lut) - hhog_mapper.ScalarVisibilityOn() - hhog_actor = vtk.vtkActor() - hhog_actor.SetMapper(hhog_mapper) - hhog_actor.SetPosition(-xmid, -ymid, -zmid) - - #cut plane - plane = vtk.vtkPlane() - plane.SetOrigin(data_scal.GetCenter()) - plane.SetNormal(0.0, 0.0, 1.0) - plane_cut = vtk.vtkCutter() - plane_cut.SetInputData(data_scal) - plane_cut.SetInputArrayToProcess( - 0, 0, 0, - vtk.vtkDataObject().FIELD_ASSOCIATION_POINTS, - vtk.vtkDataSetAttributes().SCALARS - ) - plane_cut.SetCutFunction(plane) - plane_cut.GenerateCutScalarsOff() - plane_cut.SetSortByToSortByCell() - clut = vtk.vtkLookupTable() - clut.SetHueRange(0, .67) - clut.Build() - plane_mapper = vtk.vtkPolyDataMapper() - plane_mapper.SetInputConnection(plane_cut.GetOutputPort()) - plane_mapper.ScalarVisibilityOn() - plane_mapper.SetScalarRange(data_scal.GetScalarRange()) - plane_mapper.SetLookupTable(clut) - plane_actor = vtk.vtkActor() - plane_actor.SetMapper(plane_mapper) - - #clip plane - plane_clip = vtk.vtkClipDataSet() - plane_clip.SetInputData(data_scal) - plane_clip.SetInputArrayToProcess( - 0, 0, 0, - vtk.vtkDataObject().FIELD_ASSOCIATION_POINTS, - vtk.vtkDataSetAttributes().SCALARS - ) - plane_clip.SetClipFunction(plane) - plane_clip.InsideOutOn() - clip_mapper = vtk.vtkDataSetMapper() - clip_mapper.SetInputConnection(plane_cut.GetOutputPort()) - clip_mapper.ScalarVisibilityOn() - clip_mapper.SetScalarRange(data_scal.GetScalarRange()) - clip_mapper.SetLookupTable(clut) - clip_actor = vtk.vtkActor() - clip_actor.SetMapper(clip_mapper) - - # Bounding box - outline_data = vtk.vtkOutlineFilter() - outline_data.SetInputData(data_test) - outline_mapper = vtk.vtkPolyDataMapper() - outline_mapper.SetInputConnection(outline_data.GetOutputPort()) - outline = vtk.vtkActor() - outline.SetMapper(outline_mapper) - outline.GetProperty().SetColor(0, 0, 0) - outline.GetProperty().SetLineWidth(5.0) - outline.SetPosition(-xmid, -ymid, -zmid) - - # Text - txt = vtk.vtkTextActor() - txt.SetInput('T = TEMP K') - txtprop = txt.GetTextProperty() - txtprop.SetFontFamilyToArial() - txtprop.SetFontSize(36) - txtprop.SetColor(0, 0, 0) - txt.SetDisplayPosition(30, 550) - - # Reposition the camera - ren.GetActiveCamera().Azimuth(0) - ren.GetActiveCamera().Elevation(0) - ren.GetActiveCamera().ParallelProjectionOn() - ren.GetActiveCamera().SetFocalPoint(0, 0, 0) - ren.GetActiveCamera().SetViewUp(0, 1, 0) - ren.GetActiveCamera().SetParallelScale(0.55 * ymax) - _length = max(xmax - xmin, zmax - zmin) / 2 - _height = _length / 0.26795 * 1.1 - - ren.GetActiveCamera().SetPosition(0, 0, _height) - ren.AddActor(vector) - - # Render scene - iren = vtk.vtkRenderWindowInteractor() - istyle = vtk.vtkInteractorStyleTrackballCamera() - iren.SetInteractorStyle(istyle) - iren.SetRenderWindow(ren_win) - iren.Initialize() - ren_win.Render() - screenshot(ren_win, temperature, field, plot_dir) - - @calcfunction def store_sknum_data(sky_num_out, sky_avg_num_out): """Store the skyrmion number and average skyrmion number as arrays""" @@ -419,16 +126,7 @@ def define(cls, spec): ) spec.outline( cls.submits, - if_(cls.check_for_plot_sk)( - cls.results_and_plot, - cls.combined_for_pd, - ), - if_(cls.check_for_plot_sk_number)( - cls.plot_skynumber, - cls.plot_skynumber_avg, - cls.plot_skynumber_number_heat_map, - cls.plot_sk_number_avg_heatmap, - ), + if_(cls.check_for_plot_sk)(cls.results_and_plot,), ) def check_for_plot_sk(self): @@ -440,15 +138,6 @@ def check_for_plot_sk(self): except BaseException: return False - def check_for_plot_sk_number(self): - """Check if the skyrmion number should be plotted""" - try: - if self.inputs.sk_number_plot.value > int(0): - return True - return False - except BaseException: - return False - def generate_inputs(self): """Generate the inputs for the calculation""" inputs = self.exposed_inputs( @@ -498,18 +187,6 @@ def results_and_plot(self): sm_list = [] sam_list = [] for _index_field, _ in enumerate(self.inputs.external_fields): - coord_array = self.ctx[f'T{_index_temperature}_B{_index_field}'].get_outgoing( - ).get_node_by_label('coord') - - # Size of system for plot - max_number_atoms = 1000000 - mom_states = self.ctx[f'T{_index_temperature}_B{_index_field}'].get_outgoing( - ).get_node_by_label('mom_states_traj') - points, number_atoms = read_atoms(coord_array, max_number_atoms) - vectors, colors = read_vectors_data(mom_states, number_atoms) - plot_dir = self.inputs.plot_dir.value - if self.inputs.plot_individual.value > int(0): - plot_realspace(points, vectors, colors, _index_temperature, _index_field, plot_dir) sk_data = self.ctx[f'T{_index_temperature}_B{_index_field}'].get_outgoing( ).get_node_by_label('sk_num_out') @@ -522,124 +199,3 @@ def results_and_plot(self): sk_num_for_plot = store_sknum_data(List(list=sky_num_out), List(list=sky_avg_num_out)) self.out('sk_num_for_plot', sk_num_for_plot) - - def combined_for_pd(self): - """Combine the collected data for the phase diagrame""" - if self.inputs.plot_individual.value > int(0) and self.inputs.plot_combine.value > int(0): - skyrmions_singleplot = [] - plot_dir = self.inputs.plot_dir.value - for _index_temperature, _ in enumerate(self.inputs.temperatures): - for _index_field, _ in enumerate(self.inputs.external_fields): - skyrmions_singleplot.append( - Image.open(f'{plot_dir}/T{_index_temperature}B{_index_field}.png').crop((200, 0, 2096, 1280)) - ) #we need to crop the white edges - phase_diagram = Image.new( - 'RGB', (648 * len(self.inputs.temperatures), 640 * len(self.inputs.external_fields)) - ) - for _index_temperature, _ in enumerate(self.inputs.temperatures): - for _index_field, _ in enumerate(self.inputs.external_fields): - phase_diagram.paste( - skyrmions_singleplot[len(self.inputs.external_fields) * _index_temperature + _index_field], - (0 + 648 * _index_temperature, 0 + 640 * _index_field) - ) - phase_diagram.save(f'{plot_dir}/PhaseDiagram.png', quality=100) - - def plot_skynumber(self): - """Plot the skyrmion number""" - plot_dir = self.inputs.plot_dir.value - _, axes = plt.subplots() - for _index_field, external_field in enumerate(self.inputs.external_fields): - sk_number_list = [] - for _index_temperature, _ in enumerate(self.inputs.temperatures): - sk_number_list.append( - self.ctx[f'T{_index_temperature}_B{_index_field}'].get_outgoing().get_node_by_label('sk_num_out'). - get_array('sk_number_data')[0] - ) - axes.plot(self.inputs.temperatures.get_list(), sk_number_list, label=f'B: {external_field[-6:]} T') - axes.legend(fontsize='xx-small') - axes.set_ylabel('Skyrmion number') - axes.set_xlabel('Temperature') - plt.savefig(f'{plot_dir}/sk_number.png') - - def plot_skynumber_avg(self): - """Plot the average skyrmion number""" - plot_dir = self.inputs.plot_dir.value - _, axes = plt.subplots() - for _index_field, external_field in enumerate(self.inputs.external_fields): - sk_number_list = [] - for _index_temperature, _ in enumerate(self.inputs.temperatures): - sk_number_list.append( - self.ctx[f'T{_index_temperature}_B{_index_field}'].get_outgoing().get_node_by_label('sk_num_out'). - get_array('sk_number_data')[1] - ) - axes.plot( - self.inputs.temperatures.get_list(), sk_number_list, label=f'B: {external_field[-6:]} T' - ) #this should be changed with different B - axes.legend(fontsize='xx-small') - axes.set_ylabel('Skyrmion number') - axes.set_xlabel('Temperature') - plt.savefig(f'{plot_dir}/sk_number_avg.png') - - def plot_skynumber_number_heat_map(self): - """Plot the heat map the skyrmion number""" - plot_dir = self.inputs.plot_dir.value - heat_map = [] - for _index_field, _ in enumerate(self.inputs.external_fields): - sk_number_list = [] - for _index_temperature, _ in enumerate(self.inputs.temperatures): - sk_number_list.append( - self.ctx[f'T{_index_temperature}_B{_index_field}'].get_outgoing().get_node_by_label('sk_num_out'). - get_array('sk_number_data')[0] - ) - heat_map.append(sk_number_list) - y_orginal = self.inputs.external_fields.get_list() - y_for_plot = [] - for i in y_orginal: - y_for_plot.append(float(i.split()[-1])) - - fig = go.Figure( - data=go.Heatmap( - z=heat_map, - x=self.inputs.temperatures.get_list(), - y=y_for_plot, - zsmooth='best', - zmid=0, - ) - ) - - fig.update_xaxes(range=[self.inputs.temperatures.get_list()[0], self.inputs.temperatures.get_list()[-1]]) - fig.update_yaxes(range=[y_for_plot[0], y_for_plot[-1]]) - fig.update_traces(colorscale='Jet', selector=dict(type='heatmap')) - fig.write_image(f'{plot_dir}/sk_number_heatmap.png') - - def plot_sk_number_avg_heatmap(self): - """Plot the heat map of the average skyrmion number""" - plot_dir = self.inputs.plot_dir.value - heat_map = [] - for _index_field, _ in enumerate(self.inputs.external_fields): - sk_number_list = [] - for _index_temperature, _ in enumerate(self.inputs.temperatures): - sk_number_list.append( - self.ctx[f'T{_index_temperature}_B{_index_field}'].get_outgoing().get_node_by_label('sk_num_out'). - get_array('sk_number_data')[1] - ) - heat_map.append(sk_number_list) - y_orginal = self.inputs.external_fields.get_list() - y_for_plot = [] - for i in y_orginal: - y_for_plot.append(float(i.split()[-1])) - - fig = go.Figure( - data=go.Heatmap( - z=heat_map, - x=self.inputs.temperatures.get_list(), - y=y_for_plot, - zsmooth='best', - zmid=0, - ) - ) - - fig.update_xaxes(range=[self.inputs.temperatures.get_list()[0], self.inputs.temperatures.get_list()[-1]]) - fig.update_yaxes(range=[y_for_plot[0], y_for_plot[-1]]) - fig.update_traces(colorscale='Jet', selector={'type': 'heatmap'}) - fig.write_image(f'{plot_dir}/sk_number_heatmap_avg.png') diff --git a/aiida_uppasd/workflows/temperature_cell_size.py b/aiida_uppasd/workflows/temperature_cell_size.py index cbbaae8..4101772 100644 --- a/aiida_uppasd/workflows/temperature_cell_size.py +++ b/aiida_uppasd/workflows/temperature_cell_size.py @@ -6,30 +6,14 @@ """ import os -import matplotlib.pyplot as plt -import numpy as np - from aiida import load_profile, orm -from aiida.engine import ExitCode, WorkChain +from aiida.engine import WorkChain from aiida_uppasd.workflows.temperature_restart import UppASDTemperatureRestartWorkflow load_profile() -def cal_node_query(workchain_pk): - """Simple query function that helps find output dict""" - query = orm.QueryBuilder() - query.append( - orm.WorkChainNode, - filters={'id': str(workchain_pk)}, - tag='workflow_node', - ) - query.append(orm.Dict, with_incoming='workflow_node', tag='workdict') - - return query.all() - - class MCVariableCellWorkchain(WorkChain): """This demo we only need temperature list and N(cell size) @@ -59,10 +43,7 @@ def define(cls, spec): help='plot dir ', required=False, ) - spec.outline( - cls.submit_workchains, - cls.plot_result, - ) + spec.outline(cls.submit_workchains,) def submit_workchains(self): """submit UppASDTemperatureRestartWorkflow with input""" @@ -127,55 +108,3 @@ def submit_workchains(self): future = self.submit(UppASDTemperatureRestartWorkflow, **input_uppasd) key = f'workchain_with_N_{cell_size}' self.to_context(**{key: future}) - - def plot_result(self): - """ - Basic plot function - Note that we use normalized axis like T/Tc in all figures - """ - plot_path = self.inputs.plot_dir.value - ncell_list = self.inputs.N_list - _ = self.inputs.temp_list - curie_temperature = 1043 # BCC Fe Curie temperature - - plt.figure() - _, axes = plt.subplots() - for i in ncell_list: - key = f'workchain_with_N_{i}' - sub_workchain_node = self.ctx[key] - sub_workchain_pk = sub_workchain_node.pk - self.report(f'sub_workchain {sub_workchain_pk}') - data = cal_node_query(sub_workchain_pk)[0][0].get_dict() - axes.plot( - np.array(data['temperature']) / curie_temperature, - np.array(data['magnetization']) / np.array(data['magnetization'][0]), - label=f'N={i}', - ) - axes.legend() - axes.set_xlabel('T/Tc') - axes.set_ylabel('M/Mc') - # plt.title('Temperature-magnetization') - plt.savefig(f'{plot_path}/temperature-magnetization.png') - plt.close() - - plt.figure() - _, axes = plt.subplots() - for i in ncell_list: - key = f'workchain_with_N_{i}' - sub_workchain_node = self.ctx[key] - sub_workchain_pk = sub_workchain_node.pk - self.report(f'sub_workchain {sub_workchain_pk}') - data = cal_node_query(sub_workchain_pk)[0][0].get_dict() - axes.plot( - np.array(data['temperature']) / curie_temperature, - data['energy'], - label=f'N={i}', - ) - axes.legend() - axes.set_xlabel('T/Tc') - axes.set_ylabel('Energy(mRy)') - # plt.title('temperature-energy') - plt.savefig(f'{plot_path}/temperature-energy.png') - plt.close() - - return ExitCode(0) diff --git a/aiida_uppasd/workflows/thermal_dynamic.py b/aiida_uppasd/workflows/thermal_dynamic.py index de15f94..91ad1a1 100644 --- a/aiida_uppasd/workflows/thermal_dynamic.py +++ b/aiida_uppasd/workflows/thermal_dynamic.py @@ -3,15 +3,12 @@ import json from pathlib import Path -import matplotlib.pyplot as plt -import numpy as np from numpy import asarray, gradient -import plotly.graph_objects as go from scipy import integrate from scipy.interpolate import InterpolatedUnivariateSpline from aiida.common import AttributeDict -from aiida.engine import ToContext, WorkChain, calcfunction, if_ +from aiida.engine import ToContext, WorkChain, calcfunction from aiida.orm import Dict, Int, List, Str from aiida.plugins import CalculationFactory @@ -68,58 +65,6 @@ def get_temperature_data(**kwargs): return Dict(dict=outputs) -def plot_pd( - plot_dir, - heat_map, - x_label_list, - y_label_list, - plot_name, - xlabel, - ylabel, -): - #pylint: disable=too-many-arguments - """Plotting the phase diagram""" - fig = go.Figure(data=go.Heatmap( - z=heat_map, - x=x_label_list, - y=y_label_list, - zsmooth='best', - )) - fig.update_xaxes(range=[int(float(x_label_list[0])), int(float(x_label_list[-1]))]) - fig.update_yaxes(range=[int(float(y_label_list[0])), int(float(y_label_list[-1]))]) - fig.update_traces(colorscale='Jet', selector=dict(type='heatmap')) - fig.update_layout(yaxis={'title': f'{xlabel}', 'tickangle': -90}, xaxis={'title': f'{ylabel}'}) - fig.write_image(f"{plot_dir}/{plot_name.replace(' ', '_')}.png") - - -def plot_line( - x_data, - y_data, - line_name_list, - x_label, - y_label, - plot_path, - plot_name, - leg_name_list, -): - #pylint: disable=too-many-arguments - """Line plot""" - #Since we need all lines in one plot here x and y should be a dict with the line name on that - plt.figure() - _, axes = plt.subplots() - for index, _entry in enumerate(line_name_list): - axes.plot( - x_data, - y_data[_entry], - label=f'{leg_name_list[index]}', - ) - axes.legend() - axes.set_xlabel(x_label) - axes.set_ylabel(y_label) - plt.savefig(f'{plot_path}/{plot_name}.png') - plt.close() - - class ThermalDynamicWorkflow(WorkChain): # pylint: disable=too-many-public-methods """Base workchain first #Workchain to run an UppASD simulation with automated error handling and restarts.#""" @@ -302,28 +247,6 @@ def define(cls, spec): cls.load_tasks, cls.loop_temperatures, cls.results, - if_(cls.check_mag_phase_diagram)(cls.plot_mag_phase_diagram - ).elif_(cls.plot_mag_temp)(cls.plot_mag_temp).else_(cls.error_report), - if_(cls.check_sus_phase_diagram)(cls.plot_sus_phase_diagram).elif_(cls.check_susceptibility_temp - )(cls.plot_susceptibility_temp - ).else_(cls.error_report), - if_(cls.check_cev_phase_diagram)(cls.plot_cev_phase_diagram).elif_(cls.check_specific_heat_temp - )(cls.plot_specific_heat_temp - ).else_(cls.error_report), - if_(cls.check_free_e_phase_diagram)(cls.plot_free_e_phase_diagram).elif_(cls.check_free_e_temp - )(cls.plot_free_e_temp - ).else_(cls.error_report), - if_(cls.check_entropy_phase_diagram)(cls.plot_entropy_phase_diagram).elif_(cls.check_entropy_temp - )(cls.plot_entropy_temp - ).else_(cls.error_report), - if_(cls.check_energy_phase_diagram)(cls.plot_energy_phase_diagram).elif_(cls.check_energy_temp - )(cls.plot_energy_temp - ).else_(cls.error_report), - if_(cls.check_dudt_phase_diagram)(cls.plot_dudt_phase_diagram - ).elif_(cls.check_dudt_temp)(cls.plot_dudt_temp).else_(cls.error_report), - if_(cls.check_binder_cumu_phase_diagram)(cls.plot_binder_cumu_phase_diagram - ).elif_(cls.check_binder_cumulant_temp - )(cls.plot_binder_cumulant_temp).else_(cls.error_report), ) def error_report(self): @@ -340,445 +263,6 @@ def check_mag_phase_diagram(self): return False - def plot_mag_phase_diagram(self): - """Plot the magnetization phase diagram""" - y_label_list = [] - for i in self.inputs.external_fields.get_list(): - y_label_list.append(float(max(np.array(i.split())))) - - for _, cell_size in enumerate(self.inputs.cell_size): - pd_dict = self.ctx.TD_dict[f'{cell_size}'] - pd_for_plot = [] - for i in pd_dict.keys(): - pd_for_plot.append(pd_dict[i]['magnetization']) - plot_pd( - self.inputs.plot_dir.value, pd_for_plot, self.inputs.temperatures.get_list(), y_label_list, - ('M_T' + str(cell_size)), 'B', 'T' - ) - - def check_mag_temp(self): - """Check if the magnetization phase diagram should be plotted""" - if ( - self.inputs.M_T_plot.value > int(0) and len(self.inputs.temperatures.get_list()) > 1 and - len(self.inputs.external_fields.get_list()) - ) == 1: - return True - return False - - def plot_mag_temp(self): - """Plot the magnetization as a function of temperature""" - line_name_list = [] - leg_name_list = [] - line_for_plot = {} - for _, cell_size in enumerate(self.inputs.cell_size): - line_dict = self.ctx.TD_dict[f'{cell_size}'] - leg_name_list.append(('Cell:' + cell_size.replace(' ', '_'))) - index = 0 - for i in line_dict.keys(): - line_for_plot[f"{cell_size.replace(' ', '_')}_{index}"] = line_dict[i]['magnetization'] - line_name_list.append(f"{cell_size.replace(' ', '_')}_{index}") - index = index + 1 - plot_line( - self.inputs.temperatures.get_list(), - line_for_plot, - line_name_list, - 'T', - 'M', - self.inputs.plot_dir.value, - 'M_T', - leg_name_list, - ) - - #specific_heat - def check_cev_phase_diagram(self): - """Check if the specific heat should be produced""" - if ( - self.inputs.specific_heat_phase_diagram.value > int(0) and len(self.inputs.temperatures.get_list()) > 1 and - len(self.inputs.external_fields.get_list()) - ) > 1: - return True - - return False - - def plot_cev_phase_diagram(self): - """Check if the specific heat phase diagram should be produced""" - y_label_list = [] - for i in self.inputs.external_fields.get_list(): - y_label_list.append(float(max(np.array(i.split())))) - - for _, cell_size in enumerate(self.inputs.cell_size): - pd_dict = self.ctx.TD_dict[f'{cell_size}'] - pd_for_plot = [] - for i in pd_dict.keys(): - pd_for_plot.append(pd_dict[i]['specific_heat']) - plot_pd( - self.inputs.plot_dir.value, pd_for_plot, self.inputs.temperatures.get_list(), y_label_list, - ('specific_heat_T' + str(cell_size)), 'B', 'T' - ) - - def check_specific_heat_temp(self): - """Check if the specific heat vs temperature should be produced""" - if ( - self.inputs.specific_heat_T_plot.value > int(0) and len(self.inputs.temperatures.get_list()) > 1 and - len(self.inputs.external_fields.get_list()) - ) == 1: - return True - return False - - def plot_specific_heat_temp(self): - """Plot the specific head as as function of temperature""" - line_name_list = [] - leg_name_list = [] - line_for_plot = {} - for _, cell_size in enumerate(self.inputs.cell_size): - line_dict = self.ctx.TD_dict[f'{cell_size}'] - leg_name_list.append(('Cell:' + cell_size.replace(' ', '_'))) - index = 0 - for i in line_dict.keys(): - line_for_plot[f"{cell_size.replace(' ', '_')}_{index}"] = line_dict[i]['specific_heat'] - line_name_list.append(f"{cell_size.replace(' ', '_')}_{index}") - index = index + 1 - plot_line( - self.inputs.temperatures.get_list(), line_for_plot, line_name_list, 'T', 'specific_heat', - self.inputs.plot_dir.value, 'specific_heat_T', leg_name_list - ) - - #susceptibility - def check_sus_phase_diagram(self): - """Check if the magnetic susceptibility phase diagram should be produced""" - if ( - self.inputs.susceptibility_phase_diagram.value > int(0) and len(self.inputs.temperatures.get_list()) > 1 and - len(self.inputs.external_fields.get_list()) - ) > 1: - return True - return False - - def plot_sus_phase_diagram(self): - """Plot the magnetic susceptibility phase diagram""" - y_label_list = [] - for i in self.inputs.external_fields.get_list(): - y_label_list.append(float(max(np.array(i.split())))) - - for _, cell_size in enumerate(self.inputs.cell_size): - pd_dict = self.ctx.TD_dict[f'{cell_size}'] - pd_for_plot = [] - for i in pd_dict.keys(): - pd_for_plot.append(pd_dict[i]['susceptibility']) - plot_pd( - self.inputs.plot_dir.value, pd_for_plot, self.inputs.temperatures.get_list(), y_label_list, - ('susceptibility_T' + str(cell_size)), 'B', 'T' - ) - - def check_susceptibility_temp(self): - """Check if the susceptibility vs temp should be produced""" - if ( - self.inputs.susceptibility_T_plot.value > int(0) and len(self.inputs.temperatures.get_list()) > 1 and - len(self.inputs.external_fields.get_list()) - ) == 1: - return True - return False - - def plot_susceptibility_temp(self): - """Plot the susceptibility as a function of temperature""" - line_name_list = [] - leg_name_list = [] - line_for_plot = {} - for _, cell_size in enumerate(self.inputs.cell_size): - line_dict = self.ctx.TD_dict[f'{cell_size}'] - leg_name_list.append(('Cell:' + cell_size.replace(' ', '_'))) - index = 0 - for i in line_dict.keys(): - line_for_plot[f"{cell_size.replace(' ', '_')}_{index}"] = line_dict[i]['susceptibility'] - line_name_list.append(f"{cell_size.replace(' ', '_')}_{index}") - index = index + 1 - plot_line( - self.inputs.temperatures.get_list(), line_for_plot, line_name_list, 'T', 'Susceptibility', - self.inputs.plot_dir.value, 'Susceptibility_T', leg_name_list - ) - - #free_e - def check_free_e_phase_diagram(self): - """Check if the free energy diagram should be produced""" - if ( - self.inputs.free_e_phase_diagram.value > int(0) and len(self.inputs.temperatures.get_list()) > 1 and - len(self.inputs.external_fields.get_list()) - ) > 1: - return True - return False - - def plot_free_e_phase_diagram(self): - """Plot the free energy phase diagram""" - y_label_list = [] - for i in self.inputs.external_fields.get_list(): - y_label_list.append(float(max(np.array(i.split())))) - - for _, cell_size in enumerate(self.inputs.cell_size): - pd_dict = self.ctx.TD_dict[f'{cell_size}'] - pd_for_plot = [] - for i in pd_dict.keys(): - pd_for_plot.append(pd_dict[i]['free_e']) - plot_pd( - self.inputs.plot_dir.value, - pd_for_plot, - self.inputs.temperatures.get_list(), - y_label_list, - ('free_e_T' + str(cell_size)), - 'B', - 'T', - ) - - def check_free_e_temp(self): - """Check if the free energy vs temperature should be produced""" - if ( - self.inputs.free_e_T_plot.value > int(0) and len(self.inputs.temperatures.get_list()) > 1 and - len(self.inputs.external_fields.get_list()) - ) == 1: - return True - return False - - def plot_free_e_temp(self): - """Plot the free energy as a function of temperature""" - line_name_list = [] - leg_name_list = [] - line_for_plot = {} - for _, cell_size in enumerate(self.inputs.cell_size): - line_dict = self.ctx.TD_dict[f'{cell_size}'] - leg_name_list.append(('Cell:' + cell_size.replace(' ', '_'))) - index = 0 - for i in line_dict.keys(): - line_for_plot[f"{cell_size.replace(' ', '_')}_{index}"] = line_dict[i]['free_e'] - line_name_list.append(f"{cell_size.replace(' ', '_')}_{index}") - index = index + 1 - plot_line( - self.inputs.temperatures.get_list(), line_for_plot, line_name_list, 'T', 'free_e', - self.inputs.plot_dir.value, 'free_e_T', leg_name_list - ) - - #entropy - def check_entropy_phase_diagram(self): - """Check if the entropy phase diagram should be produced""" - if ( - self.inputs.entropy_phase_diagram.value > int(0) and len(self.inputs.temperatures.get_list()) > 1 and - len(self.inputs.external_fields.get_list()) - ) > 1: - return True - return False - - def plot_entropy_phase_diagram(self): - """Plot the entropy phase diagram""" - y_label_list = [] - for i in self.inputs.external_fields.get_list(): - y_label_list.append(float(max(np.array(i.split())))) - - for _, cell_size in enumerate(self.inputs.cell_size): - pd_dict = self.ctx.TD_dict[f'{cell_size}'] - pd_for_plot = [] - for i in pd_dict.keys(): - pd_for_plot.append(pd_dict[i]['entropy']) - plot_pd( - self.inputs.plot_dir.value, pd_for_plot, self.inputs.temperatures.get_list(), y_label_list, - ('entropy_T' + str(cell_size)), 'B', 'T' - ) - - def check_entropy_temp(self): - """Check if the entropy vs temperature should be produced""" - if ( - self.inputs.entropy_T_plot.value > int(0) and len(self.inputs.temperatures.get_list()) > 1 and - len(self.inputs.external_fields.get_list()) - ) == 1: - return True - return False - - def plot_entropy_temp(self): - """Plot the entropy vs temperature""" - line_name_list = [] - leg_name_list = [] - line_for_plot = {} - for _, cell_size in enumerate(self.inputs.cell_size): - line_dict = self.ctx.TD_dict[f'{cell_size}'] - leg_name_list.append(('Cell:' + cell_size.replace(' ', '_'))) - index = 0 - for i in line_dict.keys(): - line_for_plot[f"{cell_size.replace(' ', '_')}_{index}"] = line_dict[i]['entropy'] - line_name_list.append(f"{cell_size.replace(' ', '_')}_{index}") - index = index + 1 - plot_line( - self.inputs.temperatures.get_list(), line_for_plot, line_name_list, 'T', 'Entropy', - self.inputs.plot_dir.value, 'Entropy_T', leg_name_list - ) - - #energy - def check_energy_phase_diagram(self): - """Check if the energy phase diagram should be produced""" - if ( - self.inputs.energy_phase_diagram.value > int(0) and len(self.inputs.temperatures.get_list()) > 1 and - len(self.inputs.external_fields.get_list()) - ) > 1: - return True - - return False - - def plot_energy_phase_diagram(self): - """Plot the energy phase diagram""" - y_label_list = [] - for i in self.inputs.external_fields.get_list(): - y_label_list.append(float(max(np.array(i.split())))) - - for _, cell_size in enumerate(self.inputs.cell_size): - pd_dict = self.ctx.TD_dict[f'{cell_size}'] - pd_for_plot = [] - for i in pd_dict.keys(): - pd_for_plot.append(pd_dict[i]['energy']) - plot_pd( - self.inputs.plot_dir.value, pd_for_plot, self.inputs.temperatures.get_list(), y_label_list, - ('energy_T' + str(cell_size)), 'B', 'T' - ) - - def check_energy_temp(self): - """Check if the energy vs temperature should be produced""" - if ( - self.inputs.energy_T_plot.value > int(0) and len(self.inputs.temperatures.get_list()) > 1 and - len(self.inputs.external_fields.get_list()) - ) == 1: - return True - return False - - def plot_energy_temp(self): - """Plot the energy vs temperature""" - line_name_list = [] - leg_name_list = [] - line_for_plot = {} - for _, cell_size in enumerate(self.inputs.cell_size): - line_dict = self.ctx.TD_dict[f'{cell_size}'] - leg_name_list.append(('Cell:' + cell_size.replace(' ', '_'))) - index = 0 - for i in line_dict.keys(): - line_for_plot[f"{cell_size.replace(' ', '_')}_{index}"] = line_dict[i]['energy'] - line_name_list.append(f"{cell_size.replace(' ', '_')}_{index}") - index = index + 1 - plot_line( - self.inputs.temperatures.get_list(), line_for_plot, line_name_list, 'T', 'energy', - self.inputs.plot_dir.value, 'energy_T', leg_name_list - ) - - #dudt - def check_dudt_phase_diagram(self): - """Check if the variation of energy phase diagram should be produced""" - if ( - self.inputs.dudt_phase_diagram.value > int(0) and len(self.inputs.temperatures.get_list()) > 1 and - len(self.inputs.external_fields.get_list()) - ) > 1: - return True - return False - - def plot_dudt_phase_diagram(self): - """Plot the variation of the energy phase diagram""" - y_label_list = [] - for i in self.inputs.external_fields.get_list(): - y_label_list.append(float(max(np.array(i.split())))) - - for _, cell_size in enumerate(self.inputs.cell_size): - pd_dict = self.ctx.TD_dict[f'{cell_size}'] - pd_for_plot = [] - for i in pd_dict.keys(): - pd_for_plot.append(pd_dict[i]['dudt']) - plot_pd( - self.inputs.plot_dir.value, pd_for_plot, self.inputs.temperatures.get_list(), y_label_list, - ('dudt_T' + str(cell_size)), 'B', 'T' - ) - - def check_dudt_temp(self): - """Check if the variation of energy vs temperature should be produced""" - if ( - self.inputs.dudt_T_plot.value > int(0) and len(self.inputs.temperatures.get_list()) > 1 and - len(self.inputs.external_fields.get_list()) - ) == 1: - return True - return False - - def plot_dudt_temp(self): - """Plot the variation of the energy vs temperature""" - line_name_list = [] - leg_name_list = [] - line_for_plot = {} - for _, cell_size in enumerate(self.inputs.cell_size): - line_dict = self.ctx.TD_dict[f'{cell_size}'] - leg_name_list.append(('Cell:' + cell_size.replace(' ', '_'))) - index = 0 - for i in line_dict.keys(): - line_for_plot[f"{cell_size.replace(' ', '_')}_{index}"] = line_dict[i]['dudt'] - line_name_list.append(f"{cell_size.replace(' ', '_')}_{index}") - index = index + 1 - plot_line( - self.inputs.temperatures.get_list(), - line_for_plot, - line_name_list, - 'T', - 'dudt', - self.inputs.plot_dir.value, - 'dudt_T', - leg_name_list, - ) - - #binder_cumulant - def check_binder_cumu_phase_diagram(self): - """Check if the binder cumulant phase diagram should be produced""" - if ( - self.inputs.binder_cumulant_phase_diagram.value > int(0) and - len(self.inputs.temperatures.get_list()) > 1 and len(self.inputs.external_fields.get_list()) - ) > 1: - return True - return False - - def plot_binder_cumu_phase_diagram(self): - """Plot the binder cumulant phase diagram""" - y_label_list = [] - for i in self.inputs.external_fields.get_list(): - y_label_list.append(float(max(np.array(i.split())))) - - for _, cell_size in enumerate(self.inputs.cell_size): - pd_dict = self.ctx.TD_dict[f'{cell_size}'] - pd_for_plot = [] - for i in pd_dict.keys(): - pd_for_plot.append(pd_dict[i]['binder_cumulant']) - plot_pd( - self.inputs.plot_dir.value, pd_for_plot, self.inputs.temperatures.get_list(), y_label_list, - ('binder_cumulant_T' + str(cell_size)), 'B', 'T' - ) - - def check_binder_cumulant_temp(self): - """Check if the binder cumulant as a function of temperature should be produced""" - if ( - self.inputs.binder_cumulant_T_plot.value > int(0) and len(self.inputs.temperatures.get_list()) > 1 and - len(self.inputs.external_fields.get_list()) - ) == 1: - return True - return False - - def plot_binder_cumulant_temp(self): - """Plot the binder cumulant as as function of temperature""" - line_name_list = [] - leg_name_list = [] - line_for_plot = {} - for _, cell_size in enumerate(self.inputs.cell_size): - line_dict = self.ctx.TD_dict[f'{cell_size}'] - leg_name_list.append(('Cell:' + cell_size.replace(' ', '_'))) - index = 0 - for i in line_dict.keys(): - line_for_plot[f"{cell_size.replace(' ', '_')}_{index}"] = line_dict[i]['binder_cumulant'] - line_name_list.append(f"{cell_size.replace(' ', '_')}_{index}") - index = index + 1 - plot_line( - self.inputs.temperatures.get_list(), - line_for_plot, - line_name_list, - 'T', - 'binder_cumulant', - self.inputs.plot_dir.value, - 'binder_cumulant_T', - leg_name_list, - ) - def load_tasks(self): """Load predetermined parameters for the task""" fpath = str(Path(__file__).resolve().parent.parent) + '/defaults/tasks/' @@ -844,6 +328,5 @@ def results(self): str(idx_2)].get_outgoing().get_node_by_label('cumulants') temperature_output = get_temperature_data(**outputs) td_out[f'{cell_size}'][f'{external_field}'] = temperature_output.get_dict() - self.ctx.td_dict = td_out outdict = Dict(dict=td_out).store() self.out('thermal_dynamic_output', outdict) diff --git a/pyproject.toml b/pyproject.toml index 416392a..b9c0b8a 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -42,6 +42,12 @@ pre-commit = [ "pylint~=2.15.0", "sphinx-lint~=0.6" ] +plotting = [ + "vtk>=7.0", + "plotly>=5.14.1", + "matplotlib", + "scipy" +] docs = [ "aiida-core[docs]~=2.3", "sphinx-autobuild", @@ -145,6 +151,7 @@ allowlist_externals = bash commands = bash -ec 'pre-commit run --all-files || ( git diff; git status; exit 1; )' extras = pre-commit + plotting tests [flake8]