Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Julia plotting #21

Open
wants to merge 4 commits into
base: julia_plotting
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
294 changes: 294 additions & 0 deletions maestroflame/examples/Julia_plotting.ipynb

Large diffs are not rendered by default.

525 changes: 525 additions & 0 deletions maestroflame/maestroflame/julia_env/Manifest.toml

Large diffs are not rendered by default.

3 changes: 3 additions & 0 deletions maestroflame/maestroflame/julia_env/Project.toml
Original file line number Diff line number Diff line change
@@ -0,0 +1,3 @@
[deps]
Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
PyCall = "438e738f-606a-5dbb-bf0a-cddfbfd45ab0"
148 changes: 105 additions & 43 deletions maestroflame/maestroflame/plotting.py
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,10 @@
import matplotlib
import torch
import os
import sys
import julia
import time
from numpy import errstate,isneginf, isposinf

def make_movie(file_list, save_dir='./', var='enuc', movie_name="flame.mp4"):
i = 1
Expand All @@ -16,7 +20,6 @@ def make_movie(file_list, save_dir='./', var='enuc', movie_name="flame.mp4"):
os.system("rm movie_imag*")
Video("movie.mp4", embed=True)


class plotting_standard:
#class to make it easy to plot things from driver. Set up all the data
#then just call the methods for whatever plots you want.
Expand All @@ -41,10 +44,9 @@ def __init__(self, model, fields, test_loader, cost_per_epoc,
os.mkdir(output_dir)


def do_prediction_vs_solution_plot(self):
def do_prediction_vs_solution_plot(self, use_julia=False):
############ Prediction Vs Solution Plot Should fall one y=x line.


plt.figure()
#N = react_data.output_data.shape[1]
colors = matplotlib.cm.rainbow(np.linspace(0, 1, self.nnuc+1))
Expand All @@ -66,41 +68,105 @@ def do_prediction_vs_solution_plot(self):

#for batch_idx, (data, targets) in enumerate(self.test_loader):
pred = self.model(data_whole)
#print(pred.shape)
for i in range(pred.shape[0]):
if i == 0:
for j in range(pred.shape[1]):
plt.scatter(pred[i,j], targets_whole[i,j], color=colors[j], label=self.fields[j])
else:
plt.scatter(pred[i, :], targets_whole[i, :self.nnuc+1], c=colors)


if use_julia:

outer_start_time = time.time()
start_time = time.time()


# if self.output_dir[0] != '/':
# print("You must supply output_dir as an absolute path if you are going to \
# to use julia plotting!")
# sys.exit()
#
dir_path = os.path.dirname(os.path.realpath(__file__)) #+ '/maestroflame/'


np.savetxt(self.output_dir + '/targets.txt', targets_whole)
np.savetxt(self.output_dir + '/pred.txt', pred)
#np.savetxt(self.output_dir + '/labels.txt', self.fields + [self.output_dir, dir_path], fmt="%s")

from julia.api import LibJulia
api = LibJulia.load()
#api.sysimage = "sys_plots.so"
api.sysimage = "../tools/sys.so"
api.init_julia()
from julia import Plots

print("--- Internal Timing Julia")
print("%s seconds to start julia " % (time.time() - start_time))
start_time = time.time()

#TODO: can't get labels to work.
#Plots.scatter(pred, targets_whole, labels=self.fields, dpi=600)

pred_plot = pred.cpu().detach().numpy()
targets_plot = targets_whole.cpu().detach().numpy()

Plots.scatter(pred_plot, targets_plot, dpi=600, xlabel='Prediction', ylabel='Solution')
Plots.savefig(self.output_dir + "julia_prediction_vs_solution.png")

#screen out zeros when taking log. Julia log plots can't handle zeros
#so we just do it for them.

with errstate(divide='ignore'):
pred_plot = np.log10(pred_plot)
targets_plot = np.log10(targets_plot)

#screen
pred_plot_sc = pred_plot
targets_plot_sc = targets_plot

mask = np.isinf(pred_plot) + np.isinf(targets_plot) + np.isnan(pred_plot) + np.isnan(targets_plot)

pred_plot_sc[mask]=0.
targets_plot_sc[mask]=0.
print("%s seconds to do preprocessing " % (time.time() - start_time))
start_time = time.time()

Plots.scatter(pred_plot, targets_plot, dpi=600, xlabel='log10(Prediction)', ylabel='log10(Solution)')
Plots.savefig(self.output_dir + "julia_prediction_vs_solution_log.png")
print("%s seconds to do plot " % (time.time() - start_time))
print("%s seconds total for first julia plotting method " % (time.time() - outer_start_time))


outer_start_time = time.time()
start_time = time.time()

if not os.path.isdir(dir_path + "/python_to_julia_data"):
os.mkdir(dir_path + "/python_to_julia_data")
np.savetxt(dir_path + "/python_to_julia_data/"+"pred.txt", pred)
np.savetxt(dir_path + "/python_to_julia_data/"+"targets.txt", targets_whole)
str_to_julia = self.fields + [self.output_dir]

np.savetxt(dir_path + "/python_to_julia_data/"+"labels.txt", str_to_julia, fmt="%s")
print("%s seconds to save data to file " % (time.time() - start_time))
start_time = time.time()


os.system(f"julia {dir_path}"+"/pred_v_sol.jl")
print("%s seconds to run julia pred_v_sol.jl " % (time.time() - start_time))
print("%s seconds total to run julia method 2 " % (time.time() - outer_start_time))

else:
for i in range(pred.shape[1]):
plt.scatter(pred[:, i], targets_whole[:, i], color=colors[i], label=self.fields[i])
plt.plot(np.linspace(0, 1), np.linspace(0,1), '--', color='orange')
#plt.legend(yt.load(react_data.output_files[0]).field_list, colors=colors)
plt.legend(bbox_to_anchor=(1, 1))
plt.xlabel('Prediction')
plt.ylabel('Solution')
plt.savefig(self.output_dir + "/prediction_vs_solution.png", bbox_inches='tight')

plt.yscale("log")
plt.xscale("log")
plt.savefig(self.output_dir + "/prediction_vs_solution_log.png", bbox_inches='tight')

self.model.train()
# plt.figure()
# #N = react_data.output_data.shape[1]
# colors = matplotlib.cm.rainbow(np.linspace(0, 1, self.N_fields))
# #fields = [field[1] for field in yt.load(react_data.output_files[0]).field_list]
# with torch.no_grad():
# losses = []
# for batch_idx, (data, targets) in enumerate(self.test_loader):
# pred = self.model(data)
# for i in range(pred.shape[0]):
# if i == 0 and batch_idx == 0:
# for j in range(pred.shape[1]):
# plt.scatter(pred[i,j], targets[i,j], color=colors[j], label=self.fields[j])
# else:
# plt.scatter(pred[i, :], targets[i, :], c=colors)
# # if batch_idx == 100:
# # break
plt.plot(np.linspace(0, 1), np.linspace(0,1), '--', color='orange')
#plt.legend(yt.load(react_data.output_files[0]).field_list, colors=colors)
plt.legend(bbox_to_anchor=(1, 1))
plt.xlabel('Prediction')
plt.ylabel('Solution')
plt.savefig(self.output_dir + "/prediction_vs_solution.png", bbox_inches='tight')

plt.yscale("log")
plt.xscale("log")
plt.savefig(self.output_dir + "/prediction_vs_solution_log.png", bbox_inches='tight')



def do_cost_per_epoch_plot(self):
Expand Down Expand Up @@ -223,7 +289,7 @@ def do_prediction_vs_solution_plot(self):
colors = matplotlib.cm.rainbow(np.linspace(0, 1, self.nnuc+1))
#fields = [field[1] for field in yt.load(react_data.output_files[0]).field_list]
self.model.eval()

with torch.no_grad():
losses = []

Expand All @@ -239,13 +305,9 @@ def do_prediction_vs_solution_plot(self):

#for batch_idx, (data, targets) in enumerate(self.test_loader):
pred = self.model(data_whole)
#print(pred.shape)
for i in range(pred.shape[0]):
if i == 0:
for j in range(pred.shape[1]):
plt.scatter(pred[i,j], targets_whole[i,j], color=colors[j], label=self.fields[j])
else:
plt.scatter(pred[i, :self.nnuc+1], targets_whole[i, :self.nnuc+1], c=colors)

for i in range(self.nnuc+1):
plt.scatter(pred[:, i], targets_whole[:, i], color=colors[i], label=self.fields[i])

self.model.train()

Expand Down
51 changes: 51 additions & 0 deletions maestroflame/maestroflame/pred_v_sol.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,51 @@
using Pkg
using DelimitedFiles

using Base
#Pkg.activate("julia_env")
#Pkg.add("GR")



# using Plots
#
#println(pwd())


#cut off pred_v_sol.jl
source_dir = chop(Base.source_path(), tail = 13)


pred = readdlm(source_dir * "python_to_julia_data/pred.txt")
labels = readdlm(source_dir * "python_to_julia_data/labels.txt")
targets = readdlm(source_dir * "python_to_julia_data/targets.txt")

save_dir = labels[length(labels)]
#package_dir = labels[length(labels)]
labels = labels[1:length(labels)-1]

#Pkg.activate(package_dir * "/julia_env")
#ENV["GRDIR"]=""
#Pkg.build("GR")
using Plots


labels = string.(labels)
labels = reshape(labels, (1, size(labels)[1]) )


plotd = scatter(pred, targets, labels=labels, dpi=600)
savefig(plotd, save_dir * "julia_m2_prediction_vs_solution.png")
#we now cut out data below 0. Julia plotting can't plot log plots if there are negative
# or zero data i guess
#inds = (targets .> 0) .& (pred .> 0)

#we now cut out data below 0. Julia plotting can't plot log plots if there are negative
# or zero data i guess. but to preserve shape, we just put all these points at (1,1)
inds = (targets .<= 0) .| (pred .<= 0)
targets[inds] .= 1
pred[inds] .= 1


plotd = scatter(targets, pred, xaxis=:log, yaxis=:log, labels=labels, dpi=600)
savefig(plotd, save_dir * "julia_m2_log_prediction_vs_solution.png")
82 changes: 43 additions & 39 deletions maestroflame/maestroflame/train.py
Original file line number Diff line number Diff line change
Expand Up @@ -250,22 +250,6 @@ def train(self, model, optimizer, num_epochs, criterion, save_every_N=np.Inf):
self.component_losses_train.append(component_loss/batch_idx)
self.cost_per_epoc.append(sum(plotting_losses) / len(plotting_losses))

with torch.no_grad():
if epoch % save_every_N == 0 and epoch != 0:
directory = self.output_dir+'intermediate_output/epoch'+str(epoch)+'/'
os.mkdir(directory)


torch.save(model.state_dict(), directory+'my_model.pt')
np.savetxt(directory + "/cost_per_epoch.txt", self.cost_per_epoc)
np.savetxt(directory + "/component_losses_test.txt", self.component_losses_test)
np.savetxt(directory + "/component_losses_train.txt", self.component_losses_train)


plot_class = plotting_standard(model, self.fields, self.test_loader, self.cost_per_epoc, np.array(self.component_losses_test),
np.array(self.component_losses_train), self.cost_per_epoc_test, directory)

plot_class.do_all_plots()


if self.DO_PLOTTING:
Expand Down Expand Up @@ -295,9 +279,30 @@ def train(self, model, optimizer, num_epochs, criterion, save_every_N=np.Inf):

model.train()


with torch.no_grad():
if epoch % save_every_N == 0 and epoch != 0:
directory = self.output_dir+'intermediate_output/epoch'+str(epoch)+'/'
os.mkdir(directory)


torch.save(model.state_dict(), directory+'my_model.pt')
np.savetxt(directory + "/cost_per_epoch.txt", self.cost_per_epoc)
np.savetxt(directory + "/component_losses_test.txt", self.component_losses_test)
np.savetxt(directory + "/component_losses_train.txt", self.component_losses_train)


plot_class = plotting_standard(model, self.fields, self.test_loader, self.cost_per_epoc, np.array(self.component_losses_test),
np.array(self.component_losses_train), self.cost_per_epoc_test, directory)

plot_class.do_all_plots()


self.component_losses_test = np.array(self.component_losses_test)
self.component_losses_train = np.array(self.component_losses_train)



self.model = model

if self.SAVE_MODEL:
Expand Down Expand Up @@ -642,29 +647,6 @@ def criterion_plotting(prediction, targets):
diff_losses = np.array(diff_losses)
self.different_loss_metrics.append(diff_losses.sum(axis=0)/len(losses))

with torch.no_grad():
if epoch % save_every_N == 0 and epoch != 0:
directory = self.output_dir+'intermediate_output/epoch'+str(epoch)+'/'
os.mkdir(directory)


torch.save(model.state_dict(), directory+'my_model.pt')
np.savetxt(directory + "/cost_per_epoch.txt", self.cost_per_epoc)
np.savetxt(directory + "/component_losses_test.txt", self.component_losses_test)
np.savetxt(directory + "/component_losses_train.txt", self.component_losses_train)
np.savetxt(self.output_dir + "/d_component_losses_test.txt", self.d_component_losses_test)
np.savetxt(self.output_dir + "/d_component_losses_train.txt", self.d_component_losses_train)


plot_class = plotting_pinn(self.model, self.fields, self.test_loader, self.cost_per_epoc,
np.array(self.component_losses_test), np.array(self.component_losses_train),
np.array(self.d_component_losses_test), np.array(self.d_component_losses_train),
self.cost_per_epoc_test, np.array(self.different_loss_metrics),
self.output_dir)

plot_class.do_all_plots()




if self.DO_PLOTTING:
Expand Down Expand Up @@ -715,6 +697,28 @@ def criterion_plotting(prediction, targets):
self.component_losses_test.append(component_loss/batch_idx)
self.d_component_losses_test.append(d_component_loss/batch_idx)

with torch.no_grad():
if epoch % save_every_N == 0 and epoch != 0:
directory = self.output_dir+'intermediate_output/epoch'+str(epoch)+'/'
os.mkdir(directory)


torch.save(model.state_dict(), directory+'my_model.pt')
np.savetxt(directory + "/cost_per_epoch.txt", self.cost_per_epoc)
np.savetxt(directory + "/component_losses_test.txt", self.component_losses_test)
np.savetxt(directory + "/component_losses_train.txt", self.component_losses_train)
np.savetxt(self.output_dir + "/d_component_losses_test.txt", self.d_component_losses_test)
np.savetxt(self.output_dir + "/d_component_losses_train.txt", self.d_component_losses_train)


plot_class = plotting_pinn(self.model, self.fields, self.test_loader, self.cost_per_epoc,
np.array(self.component_losses_test), np.array(self.component_losses_train),
np.array(self.d_component_losses_test), np.array(self.d_component_losses_train),
self.cost_per_epoc_test, np.array(self.different_loss_metrics),
self.output_dir)

plot_class.do_all_plots()


self.component_losses_train = np.array(self.component_losses_train)
self.component_losses_test = np.array(self.component_losses_test)
Expand Down
2 changes: 2 additions & 0 deletions maestroflame/setup.py
Original file line number Diff line number Diff line change
Expand Up @@ -9,6 +9,8 @@
license='LICENSE.txt',
description='Machine Learning interface for learning a MAESTROeX flame.',
long_description=open('README.md').read(),
package_data={'': ['pred_v_sol.jl', 'julia_env/Manifest.toml', 'julia_env/Project.toml']},
include_package_data=True,
install_requires=[
"matplotlib==3.4.1",
"optuna==2.8.0",
Expand Down
7 changes: 7 additions & 0 deletions maestroflame/tools/make_sysimage.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,7 @@
import Pkg;
Pkg.add("PackageCompiler")
Pkg.add("Plots")
Pkg.add("PyCall")
using PackageCompiler, Plots, PyCall
create_sysimage(:Plots, sysimage_path="sys_plots.so")

Loading