diff --git a/kalman/compute_errors_kalman.py b/kalman/compute_errors_kalman.py new file mode 100644 index 0000000..22a75c3 --- /dev/null +++ b/kalman/compute_errors_kalman.py @@ -0,0 +1,51 @@ +import numpy as np + +vref = [10, 250, 250, 250, 30] +#case 2 +v_no = [9.37, 246.18, 251.68, 252.37, 29.79] +v_70 = [18.21, 217.19, 198, 227.76, 29] +vnew_70 = [9.42, 244.59, 247.93, 246.89, 29.79] +v_50 = [64.05, 205.9, 195.55, 246.47, 22.58] +vnew_50 = [ 9.46 ,244.47 ,248.09 ,244.69 ,29.76] + +#acse 1 +v_no_case1 = [ 9.58 , 246.97, 252.41 , 252.47 , 29.81] +v_70_case1 = [17.93, 215.25, 193.09, 226.87, 29.02] +v_50_case1 = [64.59, 207.49, 193.08, 248.51, 22.63] +vnew_70_case1 = [ 9.51 , 237.32, 238.26, 234.33, 29.79] +vnew_50_case1 = [ 9.48, 237.58, 239.46, 231.55, 29.76] + + +def compute_relative_error(vcur,vref): + eps = 0 + for i in range(len(vref)): + eps += ((vcur[i] - vref[i])/vref[i])**2 + + return np.sqrt(eps) + + +eps_ref = compute_relative_error(v_no,vref) +eps_70 = compute_relative_error(v_70,vref) +eps_new70 = compute_relative_error(vnew_70,vref) +eps_50 = compute_relative_error(v_50,vref) +eps_new50 = compute_relative_error(vnew_50,vref) + +eps_ref_case1 = compute_relative_error(v_no_case1,vref) +eps_70_case1 = compute_relative_error(v_70_case1,vref) +eps_new70_case1 = compute_relative_error(vnew_70_case1,vref) +eps_50_case1 = compute_relative_error(v_50_case1,vref) +eps_new50_case1 = compute_relative_error(vnew_50_case1,vref) + +print('reference error = {e}'.format(e=np.round(100*eps_ref,2))) +print('venc 70% error = {e}'.format(e=np.round(100*eps_70,2))) +print('odv venc 70% error = {e}'.format(e=np.round(100*eps_new70,2))) +print('venc 50% error = {e}'.format(e=np.round(100*eps_50,2))) +print('odv venc 50% error = {e}'.format(e=np.round(100*eps_new50,2))) + +print('-------------------') + +print('reference error = {e}'.format(e=np.round(100*eps_ref_case1,2))) +print('venc 70% error = {e}'.format(e=np.round(100*eps_70_case1,2))) +print('odv venc 70% error = {e}'.format(e=np.round(100*eps_new70_case1,2))) +print('venc 50% error = {e}'.format(e=np.round(100*eps_50_case1,2))) +print('odv venc 50% error = {e}'.format(e=np.round(100*eps_new50_case1,2))) \ No newline at end of file diff --git a/kalman/graphics/figure2.py b/kalman/graphics/figure2.py index 800a762..d7881bf 100644 --- a/kalman/graphics/figure2.py +++ b/kalman/graphics/figure2.py @@ -43,7 +43,7 @@ def plot_parameters(dat, deparameterize=False, ref=None): plt.ion() - inputfile_path = 'results/aorta_leo/input.yaml' + inputfile_path = 'results/aorta/input.yaml' with open(inputfile_path) as file: inputfile = yaml.full_load(file) diff --git a/kalman/aorta.yaml b/kalman/input_files/aorta.yaml similarity index 94% rename from kalman/aorta.yaml rename to kalman/input_files/aorta.yaml index 4fcad16..f383828 100755 --- a/kalman/aorta.yaml +++ b/kalman/input_files/aorta.yaml @@ -134,18 +134,20 @@ estimation: initial_stddev: 1 - measurements: - - mesh: './meshes/coaortaH1.h5' + mesh: './meshes/coaortaH3_leo2.0.h5' fe_degree: 1 xdmf_file: 'results/aorta/measurements/u_all.xdmf' file_root: 'results/aorta/measurements/u{i}.h5' indices: 0 # indices of checkpoints to be processed. 0 == all velocity_direction: ~ - noise_stddev: 5 # standard deviation of Gaussian noise + noise_stddev: 15 # standard deviation of Gaussian noise roukf: particles: 'simplex' # unique or simplex observation_operator: 'postprocessing' #state or postprocessing reparameterize: True + ODV_functional: + enable: False + VENC: 44 # cm/s ref: 61 for 70% and 44 for 50% diff --git a/kalman/aorta_C.yaml b/kalman/input_files/aorta_C.yaml similarity index 100% rename from kalman/aorta_C.yaml rename to kalman/input_files/aorta_C.yaml diff --git a/kalman/aorta_leo.yaml b/kalman/input_files/aorta_leo.yaml similarity index 94% rename from kalman/aorta_leo.yaml rename to kalman/input_files/aorta_leo.yaml index 05eef3b..53fe255 100755 --- a/kalman/aorta_leo.yaml +++ b/kalman/input_files/aorta_leo.yaml @@ -30,22 +30,22 @@ boundary_conditions: - id: 3 type: 'windkessel' - value: [1,0,0] + value: [10,0,0] p0: [0,1333.223874] - id: 4 type: 'windkessel' - value: [100,0,0] + value: [250,0,0] p0: [0,1333.223874] - id: 5 type: 'windkessel' - value: [100,0,0] + value: [250,0,0] p0: [0,1333.223874] - id: 6 type: 'windkessel' - value: [100,0,0] + value: [250,0,0] p0: [0,1333.223874] timemarching: @@ -138,8 +138,8 @@ estimation: - mesh: './meshes/coaortaH3_leo2.0.h5' fe_degree: 1 - xdmf_file: 'results/aorta/measurements/u_all.xdmf' - file_root: 'results/aorta/measurements/u{i}.h5' + xdmf_file: 'results/aorta_leo/measurements/u_all.xdmf' + file_root: 'results/aorta_leo/measurements/u{i}.h5' indices: 0 # indices of checkpoints to be processed. 0 == all velocity_direction: ~ noise_stddev: 5 # standard deviation of Gaussian noise diff --git a/kalman/aorta_neumann.yaml b/kalman/input_files/aorta_neumann.yaml similarity index 100% rename from kalman/aorta_neumann.yaml rename to kalman/input_files/aorta_neumann.yaml diff --git a/kalman/channel2d.yaml b/kalman/input_files/channel2d.yaml similarity index 100% rename from kalman/channel2d.yaml rename to kalman/input_files/channel2d.yaml diff --git a/kalman/channel3d.yaml b/kalman/input_files/channel3d.yaml similarity index 100% rename from kalman/channel3d.yaml rename to kalman/input_files/channel3d.yaml