Skip to content
Snippets Groups Projects
Commit 10aa9858 authored by Marc Karnat's avatar Marc Karnat
Browse files

Major revamp. Reduce and unified fitting code.

parent 01ca53ec
No related branches found
No related tags found
No related merge requests found
Expand all Hide whitespace changes
Inline Side-by-side
FEM_frap_mat/Droplet_sim.ipynb
+ 23
− 31
View file @ 10aa9858
%% Cell type:code id: tags:
``` python
from math import *
from single_droplet import (solve_comparison_python_matlab, solve_python,relativeerror,
plot_relative_error_and_comp_time, plot_RMS_and_comp_time,
plot_save_result, video)
from fem_frap_mat import (solve_python,relativeerror,
plot_save_result, video,plot_error_and_comp_time)
##### PARAMETERS #####
nu = 10**-6
Xi = 7./3
a = -1. # - radius of the droplet => use 1
b = sqrt(Xi/(Xi -2))*pow(nu,1/3) # => width of Phitot and Gamma0
Partitioning = 20.
Phi_in = 0.8
gamma_in = 10
gamma_out = 230
gamma_in = 5
gamma_out = 1
u_g0 = gamma_in
e_g0 = gamma_out - gamma_in
Phi_out = Phi_in/Partitioning
e = (Phi_in - Phi_out)/2
u0 = (Phi_in + Phi_out)/2
t_final = 0.5
dt = 0.1 # Must be the same as the matlab solver if we want to compare to it
Nb_step = int(t_final/dt)
parameters = [a, b, e, u0, e_g0, u_g0, Nb_step, dt]
Plot = 1 # 1 if we want to plot, 0 otherwise
Save = 1 # 1 if we want to save the plot, 0 otherwise
Save_paraview = 1 # 1 if we want to save on paraview (.xdmf), 0 otherwise
Matlab = 1 # 1 if we want to plot Matlab Result in 'plot_save_result' , 0 otherwiseS
IC = 0 # 1: Half FRAP, 0: Full FRAP
dimension = 3 # Choose the dimension of the problem, 1: like 1D, 2 = Disk, 3 = Sphere
# This function compare the Matlab results and the Python ones:
Matlab = 1 # 1 if we want to compare with Matlab ,0 otherwise
IC = 1 # 1: Half FRAP, 0: Full FRAP
# Choose the dimension of the problem, 1: like 1D, 2 = Disk, 3 = Sphere
dimension = 3
# solve_comparison_python_matlab(Plot,Save,Save_paraview,dimension,parameters,'mesh_file.xml',
# 'save_results_file.txt','save_computing_time_file.txt','save_RMS_file.txt','save_paraview_file.xdmf')
## This function solve the python FEM:
# This function compute the Python results:
# solve_python(IC, Matlab, Plot, Save,Save_paraview, dimension, parameters,
# 'mesh.xml','save_results_file.txt','save_computing_time_file.txt',
# 'save_RMS_file.txt','save_paraview_file.xdmf')
# solve_python(IC, Plot, Save, Save_paraview, dimension, parameters,
# 'mesh_file.xml', 'save_results_file.txt','save_computing_time_file.txt','save_paraview_file.xdmf')
# Plot the result from the save file containing the python results along a line:
## Plot the result from the save file containing the python results along a line,
## we can plot or not Matlab ones:
# plot_save_result(IC,Matlab,'save_results_file1.txt','save_results_file2.txt')
# Make a video from the result from the save file containing the python results along a line:
# video(IC, 'video_save.mp4','save_results_file1.txt','save_results_file2.txt')
## Make a video from the result from the save file with the python results along a line:
# video(IC, 'video_save.mp4','save_results_file1.txt','save_results_file2.txt')
# Compute the relative error between the matlab reference and the simulation:
# relativeerror('save_results_file.txt', 'save_relative_error.txt')
## Compute the relative error between the matlab reference and the simulation:
# Plot the relative error in function of time and the mean of this plot with the computing time:
# relativeerror('save_results_file.txt', 'save_relative_error.txt')
# plot_relative_error_and_comp_time( 'save_relative_error1.txt','save_computing_time_file1.txt',
# 'save_relative_error2.txt','save_computing_time_file2.txt')
# Plot the Root Mean square results in function of time and the mean of this plot with the computing time:
## Plot the relative error or the Root Mean square resultsin function of time
## and the mean of this plot with the computing time:
## relat_or_RMS = 1 # ,=1: plot relative error, =2 plot RMS, the right file need to be put in entry
# plot_RMS_and_comp_time('save_RMS_file1.txt','save_computing_time_file1.txt','save_RMS_file2.txt',
# 'save_computing_time_file2.txt')
# plot_error_and_comp_time(relat_or_RMS,'save_relative_error.txt','save_computing_time_file.txt' )
```
%% Cell type:code id: tags:
``` python
```
......
FEM_frap_mat/droplet_fit.ipynb 0 → 100644
+ 198
− 0
View file @ 10aa9858
%% Cell type:code id: tags:
``` python
from math import *
from fit_mobility import (data_to_mesh,Artificial_data_set,error_real_data,
error_artificial_data,IC_once,fit_artificial_frap,
fit_data_from_exp, plot_conv, vector,compute_sigma,
plot_save_result,video)
r = 10.252 # radius of the droplet extract from the data
Phi_max_sim = 0.8 # = Maximum value of Phi iniside the simulation = Phi_in
background = 100 # extract from the data, use in the normalization
frame_before_FRAP = 9 # Frame before FRAP, used to initialize the simulation
frame_ic = 30 # Frame that will be used as initial conditions
frame_final = 60 # Last frame of the data used
nu = 10**-9
Xi = 7./3
a = -r # - radius of the droplet
b = sqrt(Xi/(Xi -2))*pow(nu,1/3) # width of Phitot and Gamma0
dt = 0.1 # time step of the simulation
tf = 3 # final time of the simulation
gamma_in = 100
gamma_out = 200
u_g0 = gamma_in
e_g0 = gamma_out - gamma_in
Partitioning = 20.
Phi_in = Phi_max_sim
Phi_out = Phi_in/Partitioning
e = (Phi_in - Phi_out)/2
u0 = (Phi_in + Phi_out)/2
eval_radius = 0.9 # we evaluate 90% of the domain for computing the error
kept_radius = 0.9 # we kept 90% of the domain to filter the abberations
mesh_file = 'mesh.xml'
Tif_file = 'Analysis/halfFrapBig_t000000_cropped.tif'
Mid_point_radius_file = 'Analysis/mid_point_and_radius.csv'
# file where points values used in the fitting are stocked
Point_evaluate = 'Points/Point_evaluate.txt'
# file where points X coord. used in the fitting are stocked
X_coordonate = 'Points/X_coordonate.txt'
# file where points Y coord. used in the fitting are stocked
Y_coordonate = 'Points/Y_coordonate.txt'
# file where points Z coord, used in the fitting are stocked
Z_coordonate = 'Points/Z_coordonate.txt'
save_file = 'results.txt' # file where the result along a line is saved
xdmf_file = 'test.xdmf' # file where the simulation is saved
save_parameters_file = 'save_para.txt'
plot_vector_, R_plot = vector(0,1,0,eval_radius*r) #vector(z,y,x)
save_vector_, R_save = vector(0,1,0,eval_radius*r) #vector(z,y,x)
Plot = 1 # 1 if we want to plot, 0 otherwise
Save = 0 # 1 if we want to save the plot, 0 otherwise
Save_Paraview = 0 # 1 if we want to save on paraview (.xdmf), 0 otherwise
Save_Parameters = 1 # 1 if we want to save parameters, 0 otherwise
type_data = 0
## type_data=0: initial condition extract from real data
## type_data=1: initial condition full FRAP
## type_data=2: initial condition half FRAP
noise = 0
## noise =1 : Add statistical fluctuations inside the droplet
Normalization_para = [Phi_max_sim, background]
sigma = compute_sigma(Tif_file,Normalization_para,
Mid_point_radius_file,frame_before_FRAP)
parameters = [a, b, e, u0, e_g0, u_g0, tf ,dt,Partitioning,sigma]
base_folder = 'TEST/'
x0 = [e_g0, u_g0, Partitioning]
gamma_in_seed = 1
gamma_out_seed = 1
u_g0_seed = gamma_in_seed
e_g0_seed = gamma_out_seed - gamma_in_seed
seed_parameters = [e_g0_seed,u_g0_seed,Partitioning]
## Create ref points from the data
# data_to_mesh(Plot, Save,Save_Parameters, save_parameters_file,
# Normalization_para,frame_ic,frame_final,
# frame_before_FRAP,Tif_file, Mid_point_radius_file,save_file,
# eval_radius, plot_vector_, R_plot,save_vector_, R_save,
# Point_evaluate,X_coordonate, Y_coordonate, Z_coordonate)
## Create ref points with artificial parameters,
## type_data allow to choose the initial conditions
# Artificial_data_set(type_data, noise, Plot, Save, Save_Paraview,
# Save_Parameters, save_parameters_file,parameters,
# Normalization_para,Tif_file,Mid_point_radius_file,
# frame_before_FRAP,frame_ic,eval_radius, kept_radius,
# mesh_file, save_file,xdmf_file, plot_vector_, R_plot,
# save_vector_, R_save,Point_evaluate, X_coordonate,
# Y_coordonate, Z_coordonate)
## allow to compute only one time the interpolation between
## the real data and the simulation
# Phi_u01, Phi_u1, V, v, Phi_tot_int, Phi_tot_ext = IC_once(parameters,
# Normalization_para,mesh_file,Tif_file,
# Mid_point_radius_file,frame_before_FRAP,
# frame_ic, kept_radius)
## Compute error from artificial data (Ic full FRAP or Half FRAP)
# error_artificial_data(x0,2, parameters, mesh_file, X_coordonate,Y_coordonate,
# Z_coordonate,Point_evaluate, base_folder )
## Compute error from real data or artificial data initialize with real data
# error_real_data(x0,1, parameters,Phi_tot_int, Phi_tot_ext,Phi_u01, Phi_u1,v,
# V,X_coordonate,Y_coordonate,Z_coordonate,Point_evaluate, base_folder)
## Fit mobility from artificial data (Ic full FRAP or Half FRAP)
## type_data can't be 0 here
# fit_artificial_frap(type_data,seed_parameters,parameters,mesh_file,
# X_coordonate,Y_coordonate,Z_coordonate,Point_evaluate, base_folder)
## Fit mobility from real data or artificial data initialize with real data
# fit_data_from_exp(0,seed_parameters,parameters,frame_ic, kept_radius,
# Phi_tot_int, Phi_tot_ext, Phi_u01, Phi_u1,v,V,X_coordonate,
# Y_coordonate,Z_coordonate,Point_evaluate, base_folder)
## Plot the convergence of the fitting which is save inside 'base_foler'
# plot_conv(base_folder)
## Plot save results along a line in .txt, R_save need to be in adequation
## in order to have the right plot
# plot_save_result(R_save, 'results.txt')
## Make a video of the results along a line in .txt, R_save need to be
## in adequation in order to have the right plot
# video(R_save, 'test.mp4', 'results.txt')
```
%% Cell type:code id: tags:
``` python
```
FEM_frap_mat/fem_frap_mat.py
+ 229
− 362
View file @ 10aa9858
This diff is collapsed.
FEM_frap_mat/fit_mobility.py 0 → 100644
+ 924
− 0
View file @ 10aa9858
This diff is collapsed.
Meshes/mesh_fit_quick.geo 0 → 100644
+ 34
− 0
View file @ 10aa9858
//+
SetFactory("OpenCASCADE");
Sphere(1) = {0, 0, 0, 100, -Pi/2, Pi/2, 2*Pi};
//+
lc = 20;
Point(3) = {0, 0, 0, 0.1};
Physical Surface("surfacedomain") = {1};
//+
Physical Volume("vilumedomain") = {1};
//+
//+
Field[1] = Distance;
//+
Field[1].NodesList = {3};
Field[1].NNodesByEdge = 10;
Field[2] = Threshold;
Field[2].IField = 1;
Field[2].LcMin = 2;
Field[2].LcMax = 40;
Field[2].DistMin = 11;
Field[2].DistMax = 40;
Background Field = 2;
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment