Homec4science
Diffusion Homework_sp4e_Peruzzo_SáezUribe (master) Browse

R9490/homework4/exercise_6.py
13b25592f157master

Last Change
Raw File

exercise_6.py
View Options

#!/bin/env python3
from main import main
import numpy as np
import argparse
def readPositions(planet_name, directory, nb_steps ):
"""
This function reads a csv file with the state of the planet system and it create a matrix with the positions of a given planet at different times
:param planet_name: input the planet name with a NON capital first letter
:param directory: input the directory
:return: a numpy array of sizes nb_steps X 3
"""
m = nb_steps
n = 3
planet_position = np.zeros((m, n), dtype=np.float64)
# we need this if condition because the name of the files in the Trajectories folder is different from the name
# of the files that we get out of the code and are stored in the direcotry dumps
if directory == "Trajectories/":
last = 11
elif directory == "dumps/":
last = 10
else:
print("wrong directory for reference or solution files")
for i in range(0, m):
if directory == "Trajectories/":
filename = directory + "step-" + '{:04}'.format(i) + ".csv"
elif directory == "dumps/":
filename = directory + "step-" + '{:05}'.format(i) + ".csv"
else:
print("wrong directory for reference or solution files")
with open(filename, 'r') as f:
data = f.readlines()
for line in data:
words = line.split()
if words[last] == planet_name:
planet_position[i, 0:n] = np.array(words[0:n], dtype=np.float64)
return planet_position
def computeError(positions_sol, positions_ref):
"""
This function takes two matrices with the positions of a given planet at different times and it computes the integral error.
:param positions_sol: a numpy array of sizes 365 X 3 that contains the computed coordinates XYZ of the planets at different times
:param positions_ref: a numpy array of sizes 365 X 3 that contains the reference coordinates XYZ of the planets at different times
:return: type float64 containing the integral error of the simulation
"""
# todo: check if the numer of rows of positions_sol is equal to the one of positions_ref
error = np.sqrt(np.sum(np.power(np.subtract(positions_sol, positions_ref), 2)))
return error
def generateInput(scale,planet_name,input_filename,output_filename):
"""
This function generates an input file from a given input file but by scaling velocity of a given planet.
:param scale: parameter to scale the velocity
:param planet_name: string with the name of the planet
:param input_filename: name of the file containing the input for the particle code to be scaled
:param output_filename: name of the file where to write the scaled input for the particle code
:return: nothing
"""
out_file = open(output_filename, 'w')
separator = ' '
with open(input_filename, 'r') as f:
data = f.readlines()
for line in data:
words = line.split()
if words[10] == planet_name:
scale_velocity = np.multiply(scale, np.array(words[3:6], dtype=np.float64))
string_velocity = np.array2string(scale_velocity, separator=' ') # There is a loss of precision here
out_file.write(separator.join(words[0:3]) + string_velocity[1:-2] + " " + separator.join(words[6:11]) + '\n')
else:
out_file.write(separator.join(words) + '\n')
out_file.close()
return out_file
def launchParticles(input,nb_steps,freq):
main(nb_steps,freq,input,"planet",1)
if __name__ == "__main__":
parser = argparse.ArgumentParser(description='Compute_error')
parser.add_argument('nsteps', type=int,
help='specify the number of steps to perform')
parser.add_argument('freq', type=int,
help='specify the frequency for dumps')
parser.add_argument('planet_name', type=str,
help='specify the name of the planet')
parser.add_argument('scale', type=float,
help='scale for velocity')
parser.add_argument('input_filename', type=str,
help='input filename or the complete path to the file if it is not in the same folder as the script.')
parser.add_argument('output_filename', type=str,
help='output_filename or the complete path to the file if it is not in the same folder as the script.')
args = parser.parse_args()
# getting the parameters
planet_name = args.planet_name
scale = args.scale
nb_steps = args.nsteps
freq = args.freq
# getting the name of the file containing the input for the particle code to be scaled
input_filename = args.input_filename
# getting the name of the file where to write the scaled input for the particle code
output_filename = args.output_filename
# generate the input
input = generateInput(scale, planet_name, input_filename, output_filename)
# function that gives the error for a given scaling
# velocity factor of a given planet
# todo: runAndComputeError(scale,planet_name,input,nb_steps,freq)
c4science · Help