Commit ea5e4d57 authored by Tobias Leibner's avatar Tobias Leibner

update gdt

parent 407b9ec4
Pipeline #58770 failed with stages
in 2 minutes and 43 seconds
......@@ -6,8 +6,10 @@ else
fi
rm -f $(pwd)/gdt/boltzmann.so
rm -f $(pwd)/gdt/coordinatetransformedmn.so
rm -f $(pwd)/gdt/cellmodel.so
rm -f $(pwd)/gdt/vectors.so
ln -s $(pwd)/build/${build_type}/dune-gdt/lib/libboltzmann.so $(pwd)/gdt/boltzmann.so
ln -s $(pwd)/build/${build_type}/dune-gdt/lib/libcoordinatetransformedmn.so $(pwd)/gdt/coordinatetransformedmn.so
ln -s $(pwd)/build/${build_type}/dune-gdt/lib/libcellmodel.so $(pwd)/gdt/cellmodel.so
ln -s $(pwd)/build/${build_type}/dune-gdt/lib/libvectors.so $(pwd)/gdt/vectors.so
......@@ -9,16 +9,14 @@ from hapod.boltzmann.utility import create_and_scatter_boltzmann_parameters
def visualize_solutions(grid_size=50,dimension=2):
comm_world = MPI.COMM_WORLD
parameters = create_and_scatter_boltzmann_parameters(comm_world)
solver = wrapper.Solver(dimension,
"boltzmann_snapshot_sigma_s1_%g_s2_%g_a1_%g_a2_%g" % (parameters[0], parameters[1], parameters[2],
parameters[3]), 1, grid_size, True, False, parameters[0],
parameters[1], parameters[2], parameters[3])
prefix = "boltzmann_snapshot_sigma_s1_{:g}_s2_{:g}_a1_{:g}_a2_{:g}".format(*parameters)
solver = wrapper.Solver(dimension, prefix, 1, grid_size, True, False, *parameters, False)
solver.solve()
print(parameters, " done ")
if __name__ == "__main__":
argc = len(sys.argv)
grid_size = 20 if argc < 2 else int(sys.argv[1])
grid_size = 21 if argc < 2 else int(sys.argv[1])
dimension = 2 if argc < 3 else int(sys.argv[2])
visualize_solutions(grid_size, dimension)
import resource
import sys
from timeit import default_timer as timer
import numpy as np
from pymor.algorithms.pod import pod
from hapod.boltzmann.utility import create_and_scatter_boltzmann_parameters
from hapod.coordinatetransformedmn.utility import calculate_error, create_coordinatetransformedmn_solver
from hapod.mpi import MPIWrapper
def coordinatetransformedmn_pod(grid_size, tol, logfile=None):
# get MPI communicators
mpi = MPIWrapper()
# get boltzmann solver to create snapshots
mu = create_and_scatter_boltzmann_parameters(mpi.comm_world)
solver = create_coordinatetransformedmn_solver(grid_size, mu)
# calculate Boltzmann problem trajectory
start = timer()
result = solver.solve()
mpi.comm_world.Barrier()
elapsed_data_gen = timer() - start
num_snapshots = len(result)
# gather snapshots on rank 0
result, _, total_num_snapshots, _ = mpi.comm_world.gather_on_rank_0(result, num_snapshots, num_modes_equal=False)
svals = None
# perform a POD
elapsed_pod = 0
if mpi.rank_world == 0:
result, svals = pod(result, atol=0., rtol=0., l2_err=tol * np.sqrt(total_num_snapshots))
elapsed_pod = timer() - start
# write statistics to file
if logfile is not None and mpi.rank_world == 0:
logfile.write("After the POD, there are " + str(len(result)) + " modes of " + str(total_num_snapshots) +
" snapshots left!\n")
logfile.write("The maximum amount of memory used on rank 0 was: " +
str(resource.getrusage(resource.RUSAGE_SELF).ru_maxrss / 1000.**2) + " GB\n")
elapsed = timer() - start
logfile.write("Time elapsed: " + str(elapsed) + "\n")
return result, svals, total_num_snapshots, mu, mpi, elapsed_data_gen, elapsed_pod, solver
if __name__ == "__main__":
argc = len(sys.argv)
grid_size = 21 if argc < 2 else int(sys.argv[1])
tol = 1e-3 if argc < 3 else float(sys.argv[2])
filename = "ctmn_POD_gridsize_%d_tol_%f.log" % (grid_size, tol)
logfile = open(filename, "a")
final_modes, _, total_num_snapshots, mu, mpi, _, _, _ = coordinatetransformedmn_pod(grid_size, tol * grid_size, logfile=logfile)
final_modes, win = mpi.shared_memory_bcast_modes(final_modes)
calculate_error(final_modes, grid_size, mu, total_num_snapshots, mpi, logfile=logfile)
win.Free()
logfile.close()
if mpi.rank_world == 0:
logfile = open(filename, "r")
print("\n\n\nResults:\n")
print(logfile.read())
logfile.close()
import sys
from mpi4py import MPI
from hapod.coordinatetransformedmn import wrapper
from hapod.boltzmann.utility import create_and_scatter_boltzmann_parameters
def visualize_solutions(grid_size=21):
comm_world = MPI.COMM_WORLD
parameters = create_and_scatter_boltzmann_parameters(comm_world)
prefix = "coordinatetransformedboltzmann_snapshot_sigma_s1_{:g}_s2_{:g}_a1_{:g}_a2_{:g}".format(*parameters)
solver = wrapper.Solver(prefix, 1, grid_size, True, False, *parameters)
solver.solve()
print(parameters, " done ")
if __name__ == "__main__":
argc = len(sys.argv)
grid_size = 21 if argc < 2 else int(sys.argv[1])
visualize_solutions(grid_size)
Subproject commit 58511b72b0385abe5e8251276055955a1a935298
Subproject commit cec10b07f84d19e222357825a1841fbec7d65753
from timeit import default_timer as timer
import numpy as np
from pymor.vectorarrays.numpy import NumpyVectorSpace
from hapod.coordinatetransformedmn import wrapper
def create_coordinatetransformedmn_solver(gridsize, mu):
prefix = "coordinatetransformedmn_sigma_s_s_{}_a_{}_sigma_t_s_{}_a_{}".format(*mu)
return wrapper.Solver(prefix, 2000000, gridsize, False, False, *mu)
def calculate_trajectory_error(final_modes, grid_size, mu, with_half_steps=True):
solver = create_coordinatetransformedmn_solver(grid_size, mu)
vectors = solver.solve()
vectors_np = NumpyVectorSpace(vectors[0].dim).from_data(vectors.data)
error = np.sum((vectors_np - final_modes.lincomb(vectors_np.dot(final_modes))).l2_norm()**2)
return error
def calculate_total_projection_error(final_modes,
grid_size,
mu,
total_num_snapshots,
mpi_wrapper,
with_half_steps=True):
trajectory_error = calculate_trajectory_error(final_modes, grid_size, mu, with_half_steps)
trajectory_errors = mpi_wrapper.comm_world.gather(trajectory_error, root=0)
error = 0
if mpi_wrapper.rank_world == 0:
error = np.sqrt(np.sum(trajectory_errors) / total_num_snapshots)
return error
def calculate_error(final_modes, grid_size, mu, total_num_snapshots, mpi_wrapper, with_half_steps=True, logfile=None):
''' Calculates projection error. As we cannot store all snapshots due to memory restrictions, the
problem is solved again and the error calculated on the fly'''
start = timer()
err = calculate_total_projection_error(final_modes, grid_size, mu, total_num_snapshots, mpi_wrapper,
with_half_steps)
err = err if grid_size == 0 else err / grid_size
elapsed = timer() - start
if mpi_wrapper.rank_world == 0 and logfile is not None:
logfile.write("Time used for calculating error: " + str(elapsed) + "\n")
logfile.write("l2_error is: " + str(err) + "\n")
logfile.close()
return err
This diff is collapsed.
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment