Compare commits

...

2 Commits

Author SHA1 Message Date
Aloma Blanch
e079927b84 Fixedn the time between restars plot 2020-07-16 21:34:00 -05:00
Aloma Blanch
0fb8ccae73 Added inlet waveform + saved times that will be available for ParaView 2020-07-16 21:18:08 -05:00
2 changed files with 44 additions and 9 deletions

View File

@ -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
View File

@ -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)