Georg_Brantegger/Python-DT_Slot_3
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

code cleanup and commenting I

Browse Source
This commit is contained in:
Brantegger Georg
2022-07-25 11:51:02 +02:00
parent 9204729d0b
commit 0de946f8ac
3 changed files with 99 additions and 40 deletions
Download Patch File Download Diff File Expand all files Collapse all files

102
Ausgleichsbecken/Ausgleichsbecken_class_file.py
View File

@@ -1,3 +1,4 @@
from logging import exception
import numpy as np
#importing pressure conversion function
@@ -8,8 +9,19 @@ parent = os.path.dirname(current)
sys.path.append(parent)
from functions.pressure_conversion import pressure_conversion
def FODE_function(x, h, alpha, p, rho=1000., g=9.81):
f = x*abs(x)/h*alpha+g-p/(rho*h)
def FODE_function(x,h,A,A_a,p,rho,g):
# (FODE ... first order differential equation)
# based on the outflux formula by Andreas Malcherek
# https://www.youtube.com/watch?v=8HO2LwqOhqQ
# adapted for a pressurized pipeline into which the reservoir effuses
# and flow direction
# x ... effusion velocity
# h ... level in the reservoir
# A_a ... Outflux_Area
# A ... Reservoir_Area
# g ... gravitational acceleration
# rho ... density of the liquid in the reservoir
f = x*abs(x)/h*(A_a/A-1)+g-p/(rho*h)
return f
@@ -19,67 +31,84 @@ class Ausgleichsbecken_class:
# units are used to label graphs and print units are used to have a bearable format when using pythons print()
area_unit = r'$\mathrm{m}^2$'
area_outflux_unit = r'$\mathrm{m}^2$'
density_unit = r'$\mathrm{kg}/\mathrm{m}^3$'
flux_unit = r'$\mathrm{m}^3/\mathrm{s}$'
level_unit = 'm'
pressure_unit = 'Pa'
time_unit = 's'
velocity_unit = r'$\mathrm{m}/\mathrm{s}$'
volume_unit = r'$\mathrm{m}^3$'
area_unit_print = 'm²'
area_outflux_unit_print = 'm²'
density_unit_print = 'kg/m³'
flux_unit_print = 'm³/s'
level_unit_print = 'm'
time_unit_print = 's'
pressure_unit_print = '--' # will be set by .set_pressure() method
velocity_unit_print = 'm/s'
volume_unit_print = 'm³'
g = 9.81 # m/s²
rho = 1000 # kg/m³
g = 9.81 # m/s² gravitational acceleration
# init
def __init__(self,area,outflux_area,level_min = 0,level_max = np.inf ,timestep = 1):
def __init__(self,area,outflux_area,level_min = 0,level_max = np.inf ,timestep = 1,rho = 1000):
self.area = area # base area of the rectangular structure
self.area_outflux = outflux_area # area of the outlet towards the pipeline
self.density = rho # density of the liquid in the system
self.level_min = level_min # lowest allowed water level
self.level_max = level_max # highest allowed water level
self.timestep = timestep # timestep of the simulation
# initialize for get_info
self.level = "--"
self.influx = "--"
self.level = "--"
self.outflux = "--"
self.volume = "--"
# setter
def set_volume(self):
def update_volume(self):
# sets volume in reservoir based on self.level
self.volume = self.level*self.area
def set_initial_level(self,initial_level):
self.level = initial_level
self.set_volume()
# sets the level in the reservoir and should only be called during initialization
if self.level == '--':
self.level = initial_level
self.update_volume()
else:
raise Exception('Initial level was already set once. Use the .update_level(self,timestep) method to update level based on net flux.')
def set_influx(self,influx):
# sets influx to the reservoir in m³/s
# positive influx means that liquid flows into the reservoir
self.influx = influx
def set_outflux(self,outflux):
self.outflux = outflux
self.outflux_vel = outflux/self.area_outflux
# sets outflux to the reservoir in m³/s
# positive outflux means that liquid flows out of reservoir the reservoir
self.outflux = outflux
def set_pressure(self,pressure,pressure_unit,display_pressure_unit):
def set_pressure(self,pressure,display_pressure_unit):
# sets the static pressure present at the outlet of the reservoir
# units are used to convert and display the pressure
self.pressure = pressure
self.pressure_unit = pressure_unit
self.pressure_unit_print = display_pressure_unit
def set_steady_state(self,ss_influx,ss_level,pressure_unit,display_pressure_unit):
def set_steady_state(self,ss_influx,ss_level,display_pressure_unit):
# find the steady state (ss) condition in which the net flux is zero
# set pressure acting on the outflux so that the level stays constant
ss_outflux = ss_influx
ss_outflux_vel = ss_outflux/self.area_outflux
ss_pressure = self.rho*self.g*ss_level-ss_outflux_vel**2*self.rho/2
ss_pressure = self.density*self.g*ss_level-ss_outflux_vel**2*self.density/2
self.set_initial_level(ss_level)
self.set_influx(ss_influx)
self.set_outflux(ss_outflux)
self.set_pressure(ss_pressure,pressure_unit,display_pressure_unit)
self.set_pressure(ss_pressure,display_pressure_unit)
# getter
def get_info(self, full = False):
new_line = '\n'
@@ -101,6 +130,7 @@ class Ausgleichsbecken_class:
f"Current outflux vel = {round(self.outflux_vel,3):<10} {self.velocity_unit_print} {new_line}"
f"Current pipe pressure = {round(p,3):<10} {self.pressure_unit_print} {new_line}"
f"Simulation timestep = {self.timestep:<10} {self.time_unit_print} {new_line}"
f"Density of liquid = {self.density:<10} {self.density_unit_print} {new_line}"
f"----------------------------- {new_line}")
else:
# :<10 pads the self.value to be 10 characters wide
@@ -116,9 +146,27 @@ class Ausgleichsbecken_class:
print(print_str)
def get_current_level(self):
return self.level
def get_current_influx(self):
return self.influx
def get_current_outflux(self):
return self.outflux
def get_current_volume(self):
return self.volume
def get_current_pressure(self):
return self.pressure
# methods
def update_level(self,timestep):
# update level based on net flux and timestep by calculating the volume change in
# the timestep and the converting the new volume to a level by assuming a cuboid reservoir
net_flux = self.influx-self.outflux
delta_V = net_flux*timestep
new_level = (self.volume+delta_V)/self.area
@@ -127,21 +175,25 @@ class Ausgleichsbecken_class:
def e_RK_4(self):
# update outflux and outflux velocity based on current pipeline pressure and waterlevel in reservoir
yn = self.outflux_vel
yn = self.outflux/self.area_outflux # outflux velocity
h = self.level
dt = self.timestep
p = self.pressure
# assume constant pipeline pressure during timestep (see comments in main_programm)
# assume constant pipeline pressure during timestep
# e_RK_4 timestep is way smalle than timestep of characteristic method, so this should be a valid approx.
# (furthermore I have no idea how to approximate p_hs otherwise :/ )
p_hs = self.pressure
alpha = (self.area_outflux/self.area-1)
A_a = self.area_outflux
A = self.area
h_hs = self.update_level(dt/2)
rho = self.density
g = self.g
# explicit 4 step Runge Kutta
Y1 = yn
Y2 = yn + dt/2*FODE_function(Y1, h, alpha, self.pressure)
Y3 = yn + dt/2*FODE_function(Y2, h_hs, alpha, p_hs)
Y4 = yn + dt*FODE_function(Y3, h_hs, alpha, p_hs)
ynp1 = yn + dt/6*(FODE_function(Y1, h, alpha, p)+2*FODE_function(Y2, h_hs, alpha, p_hs)+ \
2*FODE_function(Y3, h_hs, alpha, p_hs)+ FODE_function(Y4, h, alpha, p))
Y2 = yn + dt/2*FODE_function(Y1,h,A,A_a,self.pressure,rho,g)
Y3 = yn + dt/2*FODE_function(Y2,h_hs,A,A_a,p_hs,rho,g)
Y4 = yn + dt*FODE_function(Y3,h_hs,A,A_a,p_hs,rho,g)
ynp1 = yn + dt/6*(FODE_function(Y1,h,A,A_a,p,rho,g)+2*FODE_function(Y2,h_hs,A,A_a,p_hs,rho,g)+ \
2*FODE_function(Y3,h_hs,A,A_a,p_hs,rho,g)+ FODE_function(Y4,h,A,A_a,p,rho,g))
self.outflux_vel = ynp1
self.outflux = ynp1*self.area_outflux

14
Ausgleichsbecken/Ausgleichsbecken_test.ipynb
View File

@@ -2,7 +2,7 @@
"cells": [
{
"cell_type": "code",
"execution_count": 12,
"execution_count": 4,
"metadata": {},
"outputs": [],
"source": [
@@ -21,7 +21,7 @@
},
{
"cell_type": "code",
"execution_count": 13,
"execution_count": 5,
"metadata": {},
"outputs": [],
"source": [
@@ -46,7 +46,7 @@
},
{
"cell_type": "code",
"execution_count": 14,
"execution_count": 6,
"metadata": {},
"outputs": [],
"source": [
@@ -76,7 +76,7 @@
" V.pressure = pressure_vec[i]\n",
" V.e_RK_4()\n",
" V.level = V.update_level(V.timestep)\n",
" V.set_volume()\n",
" V.update_volume()\n",
" outflux_vec[i+1] = V.outflux\n",
" level_vec[i+1] = V.level\n",
" if V.level < total_min_level:\n",
@@ -87,7 +87,7 @@
},
{
"cell_type": "code",
"execution_count": 15,
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
@@ -141,7 +141,7 @@
],
"metadata": {
"kernelspec": {
"display_name": "Python 3.8.13 ('DT_Slot_3')",
"display_name": "Python 3.8.13 ('Georg_DT_Slot3')",
"language": "python",
"name": "python3"
},
@@ -160,7 +160,7 @@
"orig_nbformat": 4,
"vscode": {
"interpreter": {
"hash": "4a28055eb8a3160fa4c7e4fca69770c4e0a1add985300856aa3fcf4ce32a2c48"
"hash": "84fb123bdc47ab647d3782661abcbe80fbb79236dd2f8adf4cef30e8755eb2cd"
}
}
},

23
functions/pressure_conversion.py
View File

@@ -32,7 +32,7 @@ def pa_to_atm(p):
def pa_to_psi(p):
return p/6894.8
def pressure_conversion(pressure, input_unit = 'bar', target_unit = 'Pa'):
def pressure_conversion(pressure, input_unit = 'bar', target_unit = 'Pa', return_unit = False):
p = pressure
if input_unit.lower() == 'bar':
p_pa = bar_to_pa(p)
@@ -50,20 +50,27 @@ def pressure_conversion(pressure, input_unit = 'bar', target_unit = 'Pa'):
raise Exception('Given input unit not recognised. \n Known units are: Pa, bar, mWs, Torr, atm, psi')
if target_unit.lower() == 'bar':
return pa_to_bar(p_pa), target_unit
return_vec = [pa_to_bar(p_pa), target_unit]
elif target_unit.lower() == 'mws':
return pa_to_mWS(p_pa), target_unit
return_vec = [pa_to_mWS(p_pa), target_unit]
elif target_unit.lower() == 'torr':
return pa_to_torr(p_pa), target_unit
return_vec = [pa_to_torr(p_pa), target_unit]
elif target_unit.lower() == 'atm':
return pa_to_atm(p_pa), target_unit
return_vec = [pa_to_atm(p_pa), target_unit]
elif target_unit.lower() =='psi':
return pa_to_psi(p_pa), target_unit
return_vec = [pa_to_psi(p_pa), target_unit]
elif target_unit.lower() == 'pa':
return p_pa, target_unit
return_vec = [p_pa, target_unit]
else:
raise Exception('Given target unit not recognised. \n Known units are: Pa, bar, mWs, Torr, atm, psi')
if return_unit == True:
# return with pressure unit
return return_vec
else:
# return without pressure unit
return return_vec[0]
# testing_pressure_conversion
if __name__ == '__main__':
p = 1
@@ -72,6 +79,6 @@ if __name__ == '__main__':
for input_unit in unit_dict:
for target_unit in unit_dict:
converted_p = pressure_conversion(p,input_unit,target_unit)
converted_p = pressure_conversion(p,input_unit,target_unit,return_unit=False)
print(input_unit,target_unit)
print(converted_p)