Compare commits
2 Commits
24711c0281
...
e079927b84
Author | SHA1 | Date | |
---|---|---|---|
|
e079927b84 | ||
|
0fb8ccae73 |
37
functions.py
37
functions.py
|
@ -9,6 +9,7 @@ import pandas as pd
|
||||||
import tkinter as tk
|
import tkinter as tk
|
||||||
from statistics import mean
|
from statistics import mean
|
||||||
from scipy.signal import find_peaks
|
from scipy.signal import find_peaks
|
||||||
|
import matplotlib.transforms as mtransforms
|
||||||
|
|
||||||
|
|
||||||
def error_plot(folder,t_step,r_criteria,save):
|
def error_plot(folder,t_step,r_criteria,save):
|
||||||
|
@ -118,6 +119,40 @@ def flow(folder,N_ts,T_cyc,dt,n_cyc):
|
||||||
plt.show()
|
plt.show()
|
||||||
return Q
|
return Q
|
||||||
|
|
||||||
|
def inlet_flow_waveform(project_folder,t_btw_rst,N_ts,dt,T_cyc,n_cyc):
|
||||||
|
x = np.loadtxt(project_folder+'/ROI-1.flow')
|
||||||
|
t = x[:,0]
|
||||||
|
Q = -x[:,1]
|
||||||
|
Nt_pts = np.linspace(t_btw_rst,N_ts,int(N_ts/t_btw_rst))
|
||||||
|
t_pts = Nt_pts*dt
|
||||||
|
# Put all the time values on a single cardiac cylce
|
||||||
|
for n in range(len(t_pts)):
|
||||||
|
tmp=divmod(t_pts[n],T_cyc)
|
||||||
|
t_pts[n]=tmp[1]
|
||||||
|
if round(tmp[1],3) == 0:
|
||||||
|
t_pts[n]=T_cyc
|
||||||
|
# Interpolate the flow rate to obtain the location of the point
|
||||||
|
Q_pts = np.interp(t_pts, t, Q)
|
||||||
|
|
||||||
|
fig, ax = plt.subplots()
|
||||||
|
ax.plot(t, Q, 'r')
|
||||||
|
ax.plot(t_pts, Q_pts, 'ob')
|
||||||
|
trans_offset = mtransforms.offset_copy(ax.transData, fig=fig,
|
||||||
|
x=-0, y=0.15, units='inches')
|
||||||
|
|
||||||
|
ax.set(xlabel='Time [s]', ylabel='Flow Rate - Q [mL/s]',
|
||||||
|
title='Inlet Flow rate Waveform - 1 cycle')
|
||||||
|
ax.set_ylim([-10, 90])
|
||||||
|
ax.spines['right'].set_visible(False)
|
||||||
|
ax.spines['top'].set_visible(False)
|
||||||
|
# Adding label to the points
|
||||||
|
|
||||||
|
time = []
|
||||||
|
for i in range(0,np.unique(np.round(t_pts,3)).shape[0]):
|
||||||
|
time.append('$t_'+str(i+1)+'$')
|
||||||
|
for x, y, t in zip(t_pts[(-n_cyc-1):], Q_pts[(-n_cyc-1):], time):
|
||||||
|
plt.text(x, y, t, transform=trans_offset, fontsize=12)
|
||||||
|
plt.show()
|
||||||
|
|
||||||
|
|
||||||
|
|
16
main.py
16
main.py
|
@ -13,7 +13,7 @@ from scipy import signal
|
||||||
import statistics
|
import statistics
|
||||||
|
|
||||||
|
|
||||||
from functions import error_plot, periodicity, pressure, flow
|
from functions import error_plot, periodicity, pressure, flow, inlet_flow_waveform
|
||||||
|
|
||||||
# Selct dir
|
# Selct dir
|
||||||
Tk().withdraw()
|
Tk().withdraw()
|
||||||
|
@ -31,25 +31,25 @@ mylines = [] # Declare an empty list named mylines.
|
||||||
with open (project_folder + '/solver.inp', 'rt') as myfile: # Open lorem.txt for reading text data.
|
with open (project_folder + '/solver.inp', 'rt') as myfile: # Open lorem.txt for reading text data.
|
||||||
for myline in myfile: # For each line, stored as myline,
|
for myline in myfile: # For each line, stored as myline,
|
||||||
mylines.append(myline) # add its contents to mylines.
|
mylines.append(myline) # add its contents to mylines.
|
||||||
# print(mylines)
|
|
||||||
|
|
||||||
|
|
||||||
# Number of Timesteps - idx 3
|
# Number of Timesteps - idx 3
|
||||||
# Idea: remove the text to extract the number, the text part will be the same no matter the simulation
|
# Idea: remove the text to extract the number, the text part will be the same no matter the simulation
|
||||||
N_ts = int(mylines[3][20:-1])
|
N_ts = int(mylines[3][20:-1])
|
||||||
|
|
||||||
# Time Step Size - idx 4
|
# Time Step Size - idx 4
|
||||||
dt = float(mylines[4][16:-1])
|
dt = float(mylines[4][16:-1])
|
||||||
|
|
||||||
# Residual criteria - idx 4
|
# Residual criteria - idx 4
|
||||||
rc = float(mylines[26][18:-1])
|
rc = float(mylines[26][18:-1])
|
||||||
|
# Imesteps between Restarts - idx 6
|
||||||
|
t_btw_rst = int(mylines[6][37:-1])
|
||||||
|
|
||||||
# Cehcking convergency and periodicity
|
# Cehcking convergency and periodicity
|
||||||
error_plot(folder,dt,rc,False)
|
error_plot(folder,dt,rc,False)
|
||||||
periodicity(project,folder,dt,T_cyc,n_cyc)
|
periodicity(project,folder,dt,T_cyc,n_cyc)
|
||||||
|
|
||||||
# Pressure
|
# Pressure - Outlets
|
||||||
(DBP,MBP,SBP,PP) = pressure(folder,N_ts,T_cyc,dt,n_cyc)
|
(DBP,MBP,SBP,PP) = pressure(folder,N_ts,T_cyc,dt,n_cyc)
|
||||||
|
|
||||||
# Flow Rate
|
# Flow Rate - Outlets
|
||||||
(Q_avg) = flow(folder,N_ts,T_cyc,dt,n_cyc)
|
(Q_avg) = flow(folder,N_ts,T_cyc,dt,n_cyc)
|
||||||
|
|
||||||
|
# Inlet Flow Rate + and t saved
|
||||||
|
inlet_flow_waveform(project_folder,t_btw_rst,N_ts,dt,T_cyc,n_cyc)
|
Loading…
Reference in New Issue
Block a user