https://github.com/SarahIaquinta/adaptative_uptake_of_NPs
Tip revision: 27b6a958c5d57f02caec72211e0b6c83664e7e08 authored by Sarah Iaquinta on 24 June 2024, 06:06:49 UTC
Merge pull request #3 from anciaux/main
Merge pull request #3 from anciaux/main
Tip revision: 27b6a95
cellular_uptake_rigid_particle.py
import argparse
import pickle
import time
from functools import lru_cache
from math import exp, pi, sqrt, tan
import time
import numpy as np
import scipy.integrate
import scipy.signal
import matplotlib.pyplot as plt
import np_uptake.figures.utils as fiu
import np_uptake.model.system_definition as sysdef
import tikzplotlib
from np_uptake.figures.utils import tikzplotlib_fix_ncols
from pathlib import Path
import seaborn as sns
class EnergyComputation:
"""A class that computes the total energy of the system
Attributes:
----------
None
Methods:
-------
compute_bending_energy_zone_3(self, f, particle, s_list_region3):
Computes the variation of bending energy in region 3
compute_bending_energy_zone_2r(self, f, particle, mechanics, wrapping, membrane):
Computes the variation of bending energy in region 2r
compute_bending_energy_zone_2l(self, f, particle, mechanics, wrapping, membrane):
Computes the variation of bending energy in region 2l
sum_adimensional_bending_energy_contributions(self, adimensional_bending_energy_2r,
adimensional_bending_energy_2l, adimensional_bending_energy_3):
Computes the variation of bending energy in regions 2r, 2l and 3
compute_adimensional_adhesion_energy(self, f, particle, mechanics, wrapping, l3):
Computes the variation of adhesion energy
compute_adimensional_tension_energy(self, f, particle, mechanics, wrapping, membrane, l3, r_2r,
r_2l):
Computes the variation of tension energy
sum_adimensional_energy_contributions(self, total_adimensional_bending_energy,
adimensional_adhesion_energy, adimensional_tension_energy):
Computes the total variation of adimensional energy
compute_total_adimensional_energy_for_a_given_wrapping_degree(self, f, particle, mechanics,
wrapping, membrane):
Computes the total variation of energy by executing the intermediate functions that have
been introduced above
compute_total_adimensional_energy_during_wrapping(self, particle, mechanics, wrapping,
membrane):
Computes the total variation of energy for all the wrapping degrees defined in the method
Wrapping.wrapping_list by executing the intermediate functions that have been introduced
above
"""
def compute_bending_energy_zone_3(self, f, particle, s_list_region3):
"""Computes the variation of bending energy in region 3
Parameters:
----------
f: float
Wrapping degree
particle: class
Parameters class which defines the geometry of the particle
s_list_region3: array
Discretization of the arclength in region 3
Returns:
----------
adimensionalized_energy: float
The variation of the bending energy in the region 3 (bending of the membrane) between
the state at wrapping f and the initial state: E(f) - E(0)
"""
dpsi3_list_region3 = particle.get_squared_dpsi_region3(f)
energy = scipy.integrate.simps(dpsi3_list_region3, s_list_region3[0:-1])
adimensionalized_energy = 0.25 * particle.effective_radius * energy
return adimensionalized_energy
def compute_bending_energy_zone_2r(self, f, particle, mechanics, wrapping, membrane):
"""Computes the variation of bending energy in region 2r
Parameters:
----------
f: float
Wrapping degree
particle: class
model.system_definition.ParticleGeometry class
mechanics: class
model.system_definition.MechanicalProperties_Adaptation class
wrapping: class
model.system_definition.Wrapping class
membrane: class
model.system_definition.MembraneGeometry class
Returns:
----------
adimensionalized_energy: float
The variation of the bending energy in the region 2r (bending of the free membrane)
between the state at wrapping f and the initial state: E(f) - E(0)
"""
a = particle.get_alpha_angle(f)
t = tan(0.25 * a)
t2 = t ** 2
b = sqrt(mechanics.sigma_bar / 2)
energy = (
-8
* b
/ particle.effective_radius
* t2
* ((1 / (t2 + exp(2 * b * membrane.l2 / particle.effective_radius))) - (1 / (t2 + 1)))
)
adimensionalized_energy = 0.25 * particle.effective_radius * energy
return adimensionalized_energy
@lru_cache(maxsize=10)
def compute_bending_energy_zone_2l(self, f, particle, mechanics, wrapping, membrane):
"""Computes the variation of bending energy in region 2l
Parameters:
----------
f: float
Wrapping degree
particle: class
model.system_definition.ParticleGeometry class
mechanics: class
model.system_definition.MechanicalProperties_Adaptation class
wrapping: class
model.system_definition.Wrapping class
membrane: class
model.system_definition.MembraneGeometry class
Returns:
----------
adimensionalized_energy: float
The variation of the bending energy in the region 2l (bending of the free membrane)
between the state at wrapping f and the initial state: E(f) - E(0)
"""
adimensionalized_energy_2r = self.compute_bending_energy_zone_2r(f, particle, mechanics, wrapping, membrane)
adimensionalized_energy_2l = adimensionalized_energy_2r
return adimensionalized_energy_2l
@lru_cache(maxsize=10)
def sum_adimensional_bending_energy_contributions(
self,
adimensional_bending_energy_2r,
adimensional_bending_energy_2l,
adimensional_bending_energy_3,
):
"""Computes the variation of bending energy in regions 2r, 2l and 3
Parameters:
----------
adimensional_bending_energy_2r: float
The variation of the bending energy in the region 2r (bending of the membrane) between
the state at wrapping f and the initial state: E(f) - E(0)
adimensional_bending_energy_2l: float
The variation of the bending energy in the region 2l (bending of the membrane) between
the state at wrapping f and the initial state: E(f) - E(0)
adimensional_bending_energy_3: float
The variation of the bending energy in the region 3 (bending of the membrane) between
the state at wrapping f and the initial state: E(f) - E(0)
Returns:
----------
total_adimensional_bending_energy: float
The variation of the bending energy in regions 2r, 2l and 3 between the state at
wrapping f and the initial state: E(f) - E(0)
"""
total_adimensional_bending_energy = (
adimensional_bending_energy_2r + adimensional_bending_energy_2l + adimensional_bending_energy_3
)
return total_adimensional_bending_energy
@lru_cache(maxsize=10)
def compute_adimensional_adhesion_energy(self, f, particle, mechanics, wrapping, l3):
"""Computes the variation of adhesion energy
Parameters:
----------
f: float
Wrapping degree
particle: class
model.system_definition.ParticleGeometry class
mechanics: class
model.system_definition.MechanicalProperties_Adaptation class
wrapping: class
model.system_definition.Wrapping class
l3: float
total arclength of region 3
Returns:
----------
adimensional_adhesion_energy: float
The variation of the adhesion energy
between the state at wrapping f and the initial state: E(f) - E(0)
"""
adimensional_adhesion_energy = -mechanics.gamma_bar(f, wrapping) * l3 * 0.25 / particle.effective_radius
return adimensional_adhesion_energy
def compute_adimensional_tension_energy(self, f, particle, mechanics, wrapping, membrane, l3, r_2r, r_2l):
"""Computes the variation of tension energy
Parameters:
----------
f: float
Wrapping degree
particle: class
model.system_definition.ParticleGeometry class
mechanics: class
model.system_definition.MechanicalProperties_Adaptation class
wrapping: class
model.system_definition.Wrapping class
membrane: class
model.system_definition.MembraneGeometry class
l3: float
total arclength of region 3
r2r: array
r abscises in region 2r
r2l: array
r abscises in region 2l
Returns:
----------
adimensional_tension_energy: float
The variation of the tension energy
between the state at wrapping f and the initial state: E(f) - E(0)
"""
adimensional_tension_energy = (
mechanics.sigma_bar * (l3 + 2 * membrane.l2 - (r_2r[-1] - r_2l[-1])) * 0.25 / particle.effective_radius
)
return adimensional_tension_energy
@lru_cache(maxsize=10)
def sum_adimensional_energy_contributions(
self, total_adimensional_bending_energy, adimensional_adhesion_energy, adimensional_tension_energy
):
"""Computes the total variation of adimensional energy
Parameters:
----------
adimensional_bending_energy: float
The variation of the bending energy
between the state at wrapping f and the initial state: E(f) - E(0)
adimensional_adhesion_energy: float
The variation of the adhesion energy
between the state at wrapping f and the initial state: E(f) - E(0)
adimensional_tension_energy: float
The variation of the tension energy
between the state at wrapping f and the initial state: E(f) - E(0)
Returns:
----------
total_adimensional_energy: float
The sum of the contributions of energy
"""
total_adimensional_energy = (
total_adimensional_bending_energy + adimensional_adhesion_energy + adimensional_tension_energy
)
return total_adimensional_energy
@lru_cache(maxsize=10)
def compute_total_adimensional_energy_for_a_given_wrapping_degree(self, f, particle, mechanics, wrapping, membrane):
"""Computes the total variation of energy by executing the intermediate functions that have
been introduced above
Parameters:
----------
f: float
Wrapping degree
particle: class
model.system_definition.ParticleGeometry class
mechanics: class
model.system_definition.MechanicalProperties_Adaptation class
wrapping: class
model.system_definition.Wrapping class
membrane: class
model.system_definition.MembraneGeometry class
Returns:
----------
total_adimensional_energy: float
The variation of the total adimensional energy
"""
_, _, _, _, s_list_region3, l3, _, _ = particle.define_particle_geometry_variables(f)
r_2r, _, r_2l, _ = membrane.compute_r2r_r2l_z2r_z2l_from_analytic_expression(f, particle, mechanics, wrapping)
adimensional_bending_energy_2r = self.compute_bending_energy_zone_2r(f, particle, mechanics, wrapping, membrane)
adimensional_bending_energy_2l = self.compute_bending_energy_zone_2l(f, particle, mechanics, wrapping, membrane)
adimensional_bending_energy_3 = self.compute_bending_energy_zone_3(f, particle, s_list_region3)
total_adimensional_bending_energy = self.sum_adimensional_bending_energy_contributions(
adimensional_bending_energy_2r, adimensional_bending_energy_2l, adimensional_bending_energy_3
)
adimensional_adhesion_energy = self.compute_adimensional_adhesion_energy(f, particle, mechanics, wrapping, l3)
adimensional_tension_energy = self.compute_adimensional_tension_energy(
f, particle, mechanics, wrapping, membrane, l3, r_2r, r_2l
)
total_adimensional_energy = self.sum_adimensional_energy_contributions(
total_adimensional_bending_energy, adimensional_adhesion_energy, adimensional_tension_energy
)
return total_adimensional_energy
def compute_total_adimensional_energy_during_wrapping(self, particle, mechanics, wrapping, membrane):
"""Computes the total variation of energy for all the wrapping degrees defined in the
method Wrapping.wrapping_list by executing the intermediate functions that have been
introduced above
Parameters:
----------
particle: class
model.system_definition.ParticleGeometry class
mechanics: class
model.system_definition.MechanicalProperties_Adaptation class
wrapping: class
model.system_definition.Wrapping class
membrane: class
model.system_definition.MembraneGeometry class
Returns:
----------
adimensional_total_energy_variation_list: array
The variation of the total adimensional energy for the wrapping degrees defined in
Wrapping.wrapping_list
energy_variation_computation_time_list: array
Computation time compute each element of adimensional_total_energy_variation_list
"""
adimensional_total_energy_variation_list = np.zeros_like(wrapping.wrapping_list)
energy_variation_computation_time_list = np.zeros_like(adimensional_total_energy_variation_list)
for i in range(len(wrapping.wrapping_list)):
f = wrapping.wrapping_list[i]
start = time.process_time()
total_adimensional_energy = self.compute_total_adimensional_energy_for_a_given_wrapping_degree(
f, particle, mechanics, wrapping, membrane
)
end = time.process_time()
energy_variation_computation_time = end - start
adimensional_total_energy_variation_list[i] = total_adimensional_energy
energy_variation_computation_time_list[i] = energy_variation_computation_time
return adimensional_total_energy_variation_list, energy_variation_computation_time_list
def determine_eq_energy(f_list, energy_list):
min_energy_index_list = scipy.signal.argrelmin(energy_list)
max_energy_index_list = scipy.signal.argrelmax(energy_list)
min_energy_index_list = list(min_energy_index_list[0])
max_energy_index_list = list(max_energy_index_list[0])
min_energy_index_list_initial = min_energy_index_list.copy()
max_energy_index_list_initial = max_energy_index_list.copy()
if len(max_energy_index_list_initial) > 0:
for i in range(min(len(min_energy_index_list_initial), len(max_energy_index_list_initial))):
index_min = min_energy_index_list_initial[i]
index_max = max_energy_index_list_initial[i]
diff_index = abs(index_min - index_max)
close_indices = diff_index == 1
if close_indices: # check if the minimum is directly after of before a minimum
# (i.e. if there is a peak of energy due to artefacts)
min_energy_index_list.remove(index_min)
max_energy_index_list.remove(index_max)
# check if the minimum is reached for f_list[-1]
if energy_list[-1] < energy_list[-2]:
min_energy_index_list = min_energy_index_list + [-1]
# check if the minimum is reached for f_list[0]
if energy_list[0] < energy_list[1]:
min_energy_index_list = [f_list[0]] + min_energy_index_list
if len(min_energy_index_list) == 0:
min_energy_index_list = [0]
min_energy_list = [energy_list[int(k)] for k in min_energy_index_list]
f_min_energy_list = [f_list[int(k)] for k in min_energy_index_list]
if len(max_energy_index_list) > 0:
max_energy_list = [energy_list[int(k)] for k in max_energy_index_list]
f_max_energy_list = [f_list[int(k)] for k in max_energy_index_list]
else:
max_energy_list = []
f_max_energy_list = []
# managing possible scipy.signal.argrelextrema outuput types
if type(min_energy_list[0]) == np.ndarray:
min_energy_list = min_energy_list[0]
f_min_energy_list = f_min_energy_list[0]
f_eq = f_min_energy_list[0]
energy_eq = min_energy_list[0]
return f_eq, energy_eq
def determine_min_energy(f_list, energy_list):
min_energy_index_list = scipy.signal.argrelmin(energy_list)
max_energy_index_list = scipy.signal.argrelmax(energy_list)
min_energy_index_list = list(min_energy_index_list[0])
max_energy_index_list = list(max_energy_index_list[0])
min_energy_index_list_initial = min_energy_index_list.copy()
max_energy_index_list_initial = max_energy_index_list.copy()
if len(max_energy_index_list_initial) > 0:
for i in range(min(len(min_energy_index_list_initial), len(max_energy_index_list_initial))):
index_min = min_energy_index_list_initial[i]
index_max = max_energy_index_list_initial[i]
diff_index = abs(index_min - index_max)
close_indices = diff_index == 1
if close_indices: # check if the minimum is directly after of before a minimum
# (i.e. if there is a peak of energy due to artefacts)
min_energy_index_list.remove(index_min)
max_energy_index_list.remove(index_max)
# check if the minimum is reached for f_list[-1]
if energy_list[-1] < energy_list[-2]:
min_energy_index_list = min_energy_index_list + [-1]
# check if the minimum is reached for f_list[0]
if energy_list[0] < energy_list[1]:
min_energy_index_list = [f_list[0]] + min_energy_index_list
if len(min_energy_index_list) == 0:
min_energy_index_list = [0]
min_energy_list = [energy_list[int(k)] for k in min_energy_index_list]
f_min_energy_list = [f_list[int(k)] for k in min_energy_index_list]
if len(max_energy_index_list) > 0:
max_energy_list = [energy_list[int(k)] for k in max_energy_index_list]
f_max_energy_list = [f_list[int(k)] for k in max_energy_index_list]
else:
max_energy_list = []
f_max_energy_list = []
# managing possible scipy.signal.argrelextrema outuput types
if type(min_energy_list[0]) == np.ndarray:
min_energy_list = min_energy_list[0]
f_min_energy_list = f_min_energy_list[0]
f_eq = f_min_energy_list[0]
energy_eq = min_energy_list[0]
return f_min_energy_list, min_energy_list
def plot_energy(
particle, mechanics, membrane, wrapping, energy_computation, createfigure, fonts, xticks, xticklabels, savefigure
):
"""Plots the evolution of the adimensional variation of energy_computation during wrapping
Parameters:
----------
particle: class
model.system_definition.ParticleGeometry class
mechanics: class
model.system_definition.MechanicalProperties_Adaptation class
membrane: class
model.system_definition.MembraneGeometry class
wrapping: class
model.system_definition.Wrapping class
createfigure: class
figures.utils.CreateFigure class
fonts: class
figures.utils.Fonts class
xticks: class
figures.utils.XTicks class
xticklabels: class
figures.utils.XTickLabels class
savefigure: class
figures.utils.SaveFigure class
Returns:
-------
None
"""
energy_list, _ = energy_computation.compute_total_adimensional_energy_during_wrapping(
particle, mechanics, wrapping, membrane
)
fig = createfigure.square_figure_7(pixels=360)
ax = fig.gca()
ax.plot(
wrapping.wrapping_list,
energy_list,
"-k",
label=r"$\overline{r} = $"
+ str(np.round(particle.r_bar, 2))
+ r" ; $\overline{\gamma}_0 = $"
+ str(mechanics.gamma_bar_0)
+ r" ; $\overline{\sigma} = $"
+ str(mechanics.sigma_bar),
linewidth=4,
)
ax.set_xticks(xticks.energy_plots())
ax.set_xticklabels(
xticklabels.energy_plots(),
font=fonts.serif(),
fontsize=fonts.axis_legend_size(),
)
ax.set_yticklabels(
ax.get_yticks(),
font=fonts.serif(),
fontsize=fonts.axis_legend_size(),
)
ax.legend(prop=fonts.serif(), loc="lower left", framealpha=0.9)
ax.set_xlabel(r"$f$ [ - ]", font=fonts.serif(), fontsize=fonts.axis_label_size())
ax.set_ylabel(r"$\overline{\Delta E}$ [ - ]", font=fonts.serif(), fontsize=fonts.axis_label_size())
savefigure.save_as_png(fig, "DeltaE_vs_f")
tikzplotlib_fix_ncols(fig)
current_path = Path.cwd()
tikzplotlib.save(current_path/"DeltaE_vs_f.tex")
def plot_energy_article(
wrapping, energy_computation, createfigure, fonts, xticks, xticklabels, savefigure
):
"""Plots the evolution of the adimensional variation of energy_computation during wrapping
Parameters:
----------
particle: class
model.system_definition.ParticleGeometry class
mechanics: class
model.system_definition.MechanicalProperties_Adaptation class
membrane: class
model.system_definition.MembraneGeometry class
wrapping: class
model.system_definition.Wrapping class
createfigure: class
figures.utils.CreateFigure class
fonts: class
figures.utils.Fonts class
xticks: class
figures.utils.XTicks class
xticklabels: class
figures.utils.XTickLabels class
savefigure: class
figures.utils.SaveFigure class
Returns:
-------
None
"""
mechanics_article = sysdef.MechanicalProperties_Adaptation(
testcase="testcase",
gamma_bar_r=1,
gamma_bar_fs=0,
gamma_bar_lambda=1,
gamma_bar_0=6,
sigma_bar=2,
)
particle_article = sysdef.ParticleGeometry(r_bar=0.3, particle_perimeter=args.particle_perimeter, sampling_points_circle=300)
membrane_article = sysdef.MembraneGeometry(particle_article, sampling_points_membrane=100)
energy_list_article, _ = energy_computation.compute_total_adimensional_energy_during_wrapping(
particle_article, mechanics_article, wrapping, membrane_article
)
f_eq, energy_eq = determine_eq_energy(wrapping.wrapping_list, energy_list_article)
f_min_energy_list, min_energy_list = determine_min_energy(wrapping.wrapping_list, energy_list_article)
f_eq, wrapping_phase_number, wrapping_phase, energy_list, time_list = identify_wrapping_phase(
particle_article, mechanics_article, membrane_article, wrapping, energy_computation
)
palette = sns.color_palette("Greens", 3)
fig = createfigure.square_figure_7(pixels=360)
ax = fig.gca()
ax.plot(
wrapping.wrapping_list,
energy_list_article,
"-k",
label=r"$\overline{r} = $"
+ str(np.round(particle.r_bar, 2))
+ r" ; $\overline{\gamma}_0 = $"
+ str(mechanics.gamma_bar_0)
+ r" ; $\overline{\sigma} = $"
+ str(mechanics.sigma_bar),
linewidth=4,
)
ax.plot(f_eq, energy_eq, marker="o", markersize=15, markeredgecolor="r", markeredgewidth=2, fillstyle="none")
for i in range(len(f_min_energy_list)):
ax.plot(f_min_energy_list[i], min_energy_list[i], marker="o", markersize=10, markeredgecolor="k", markerfacecolor=palette[0])
ax.set_xticks(xticks.energy_plots())
ax.set_xticklabels(
xticklabels.energy_plots(),
font=fonts.serif(),
fontsize=fonts.axis_legend_size(),
)
ax.set_yticks([-2, -1, 0, 1])
ax.set_yticklabels(
[-2, -1, 0, 1],
font=fonts.serif(),
fontsize=fonts.axis_legend_size(),
)
ax.set_ylim((-2.2, 1.4))
# ax.legend(prop=fonts.serif(), loc="lower left", framealpha=0.9)
ax.set_xlabel(r"$f$ [ - ]", font=fonts.serif(), fontsize=fonts.axis_label_size())
ax.set_ylabel(r"$\overline{\Delta E}$ [ - ]", font=fonts.serif(), fontsize=fonts.axis_label_size())
savefigure.save_as_png(fig, "DeltaE_vs_f_article_fig2a")
tikzplotlib_fix_ncols(fig)
current_path = Path.cwd()
tikzplotlib.save(current_path/"DeltaE_vs_f_article_fig2a.tex")
def plot_np_membrane_wrapping(f, particle, mechanics, membrane, wrapping, createfigure, fonts, savefigure):
r2r_list, z2r_list, r2l_list, z2l_list = membrane.compute_r2r_r2l_z2r_z2l_from_analytic_expression(f, particle, mechanics, wrapping)
r_list_region_1, z_list_region_1, r_list_region_3, z_list_region_3 = particle.compute_r_z_list(f)
fig = createfigure.rectangle_figure(pixels=360)
palette = sns.color_palette("Paired")
color_paired_dark_blue = palette[1]
color_paired_dark_red = palette[5]
ax = fig.gca()
ax.plot(r2r_list, z2r_list, '-', color=color_paired_dark_red, lw=2)
ax.plot(r2l_list, z2l_list, '-', color=color_paired_dark_red, lw=2)
ax.plot(r_list_region_1, z_list_region_1, '-', color=color_paired_dark_blue, lw=2)
ax.plot(r_list_region_3, z_list_region_3, '-k')
ax.set_aspect("equal", adjustable="box")
ax.set_xticks([-10, 0, 10])
ax.set_xticklabels(
[-10, 0, 10],
font=fonts.serif_3horizontal(),
fontsize=24,
)
ax.set_yticks([-2, 0, 2, 4])
ax.set_yticklabels(
[-2, 0, 2, 4],
font=fonts.serif_3horizontal(),
fontsize=24,
)
ax.set_ylim((-3, 4.1))
ax.set_xlim((-10, 10))
ax.set_xlabel("r " + r"$( \times 100)$ [nm]", font=fonts.serif(), fontsize=24)
ax.set_ylabel("z " + r"$( \times 100)$ [nm]", font=fonts.serif(), fontsize=22)
ax.legend(prop=fonts.serif(), loc='lower right', framealpha=0.7)
savefigure.save_as_png(fig, "system_during_wrapping")
def identify_wrapping_phase(particle, mechanics, membrane, wrapping, energy_computation):
"""Identifies the wrapping phase following the process introduced in [1]
Parameters:
----------
particle: class
model.system_definition.ParticleGeometry object
mechanics: class
model.system_definition.MechanicalProperties_Adaptation object
membrane: class
model.system_definition.MembraneGeometry object
energy_computation: class
model.system_definition.EnergyComputation object
Returns:
-------
f_eq: float
Wrapping degree at equilibrium
wrapping_phase_number: float
Phase number (1, 2 or 3)
wrapping_phase: str
The wrapping phase as an intelligible string
energy_list: array
The variation of the total adimensional energy for the wrapping degrees defined in
Wrapping.wrapping_list
Output from function EnergyComputation.compute_total_adimensional_energy_during_wrapping
time_list: array
Computation time compute each element of adimensional_total_energy_variation_list
Output from function EnergyComputation.compute_total_adimensional_energy_during_wrapping
"""
pickle.dumps(energy_computation.compute_total_adimensional_energy_during_wrapping)
pickle.dumps(membrane.compute_r2r_r2l_z2r_z2l_from_analytic_expression)
energy_list, time_list = energy_computation.compute_total_adimensional_energy_during_wrapping(
particle, mechanics, wrapping, membrane
)
min_energy_index_list = scipy.signal.argrelmin(energy_list)
min_energy_index_list = min_energy_index_list[0]
# check if the minimum is reached for wrapping.wrapping_list[-1]
if energy_list[-1] < energy_list[-2]:
min_energy_index_list = np.concatenate((min_energy_index_list, np.array([-1])), axis=None)
# check if the minimum is reached for wrapping.wrapping_list[0]
if energy_list[0] < energy_list[1]:
min_energy_index_list = np.concatenate((np.array(wrapping.wrapping_list[0]), min_energy_index_list), axis=None)
if len(min_energy_index_list) == 0:
min_energy_index_list = [0]
min_energy_list = [energy_list[int(k)] for k in min_energy_index_list]
f_min_energy_list = [wrapping.wrapping_list[int(k)] for k in min_energy_index_list]
# managing possible scipy.signal.argrelextrema outuput types
if type(min_energy_list[0]) == np.ndarray:
min_energy_list = min_energy_list[0]
f_min_energy_list = f_min_energy_list[0]
f_eq = f_min_energy_list[0]
wrapping_phase_number = 0
wrapping_phase = "0"
if f_eq < 0.2: # check if wrapping phase is phase 1, according to [1]
wrapping_phase_number = 1
wrapping_phase = "no wrapping"
else:
r2r, _, r2l, _ = membrane.compute_r2r_r2l_z2r_z2l_from_analytic_expression(f_eq, particle, mechanics, wrapping)
intersection_membrane = min(r2r) - max(r2l)
wrapping_phase_number = 3 if intersection_membrane < 0 else 2
wrapping_phase = "full wrapping" if intersection_membrane < 0 else "partial wrapping"
return (f_eq, wrapping_phase_number, wrapping_phase, energy_list, time_list)
def parse_arguments():
"""Parses arguments to run the code in terminal
Parameters:
----------
None
Returns:
-------
args: class
Parse of the arguments of the code
"""
parser = argparse.ArgumentParser()
parser.add_argument(
"-gamma_bar_0",
"--gamma_bar_0",
required=False,
default=10.0,
type=float,
help="initial adimensional lineic adhesion bt membrane and NP. Default value = 10.",
)
parser.add_argument(
"-gamma_bar_A",
"--gamma_bar_r",
required=False,
default=1.0,
type=float,
help="ratio of adimensional lineic adhesion bt membrane and NP. Default value = 1.",
)
parser.add_argument(
"-gamma_bar_D",
"--gamma_bar_fs",
required=False,
default=0.0,
type=float,
help="adimensional lineic adhesion between the membrane and the particle, inflexion point. Default value = 0.",
)
parser.add_argument(
"-gamma_bar_S",
"--gamma_bar_lambda",
required=False,
default=10.0,
type=float,
help="adimensional lineic adhesion between the membrane and the particle, smoothness. Default value = 10.",
)
parser.add_argument(
"-sigma_bar_0",
"--sigma_bar_0",
required=False,
default=2.0,
type=float,
help="adimensional membrane tension, value before wrapping. Default value = 2.",
)
parser.add_argument(
"-r_bar", "--r_bar", required=False, default=1.0, type=float, help="particle aspect ratio. Default value = 1."
)
parser.add_argument(
"-p",
"--particle_perimeter",
required=False,
default=2 * pi,
type=float,
help="particle perimeter. Default value = 2pi.",
)
args = parser.parse_args()
return args
if __name__ == "__main__":
start = time.time()
args = parse_arguments()
createfigure = fiu.CreateFigure()
fonts = fiu.Fonts()
savefigure = fiu.SaveFigure()
xticks = fiu.XTicks()
xticklabels = fiu.XTickLabels()
particle = sysdef.ParticleGeometry(r_bar=args.r_bar, particle_perimeter=args.particle_perimeter, sampling_points_circle=300)
mechanics = sysdef.MechanicalProperties_Adaptation(
testcase="testcase",
gamma_bar_r=args.gamma_bar_r,
gamma_bar_fs=args.gamma_bar_fs,
gamma_bar_lambda=args.gamma_bar_lambda,
gamma_bar_0=args.gamma_bar_0,
sigma_bar=args.sigma_bar_0,
)
membrane = sysdef.MembraneGeometry(particle, sampling_points_membrane=100)
wrapping = sysdef.Wrapping(wrapping_list=np.arange(0.03, 0.97, 0.003125))
energy_computation = EnergyComputation()
plot_energy(
particle,
mechanics,
membrane,
wrapping,
energy_computation,
createfigure,
fonts,
xticks,
xticklabels,
savefigure,
)
plot_energy_article(
wrapping,
energy_computation,
createfigure,
fonts,
xticks,
xticklabels,
savefigure,
)
f_eq, wrapping_phase_number, wrapping_phase, energy_list, time_list = identify_wrapping_phase(
particle, mechanics, membrane, wrapping, energy_computation
)
print("wrapping degree at equilibrium = ", np.round(f_eq, 2))
print("wrapping phase at equilibrium: ", wrapping_phase)
print('time',time.time() - start)
f = f_eq
plot_np_membrane_wrapping(f, particle, mechanics, membrane, wrapping, createfigure, fonts, savefigure)
