VesicleForms.jl

This package provides routines to solve the shape equations for axi-symmetric liquid membranes.

Functions

VesicleForms.VesicleType

Basic data structure to hold all parameters that uniquely define a vesicle.

  • spont_curvature::Real

  • area::Real

  • volume::Real

  • parameters::Array{F,1} where F<:Real

    these parameters are: osmotic pressure, surface tension, s-pole curvature, n-pole curvature, arc length

  • red_volume::Real

  • bend_rigidity::Real

source
VesicleForms.ODEshoot_root!Function

Finds the best parameters of the vesicle using root-finding.

ODEshoot_root!(vesicle, tspan)
ODEshoot_root!(vesicle, tspan, residual; solver, kwargs...)

Positional arguments

  • vesicle : input vesicle (gets mutated)
  • tspan : the interval of the arclength to integrate over, typically (0.0, 1.0)
  • residual: function that calculates the residual for the root-finding algorithm (root_residualAV)

Keyword arguments

  • solver : Currently supported solvers are :nlsolve and :minpack (:nlsolve)
  • additional keywords are passed to the solver routine
source
VesicleForms.ODEshoot_root!Function
ODEshoot_root!(, vesicle, tspan)
ODEshoot_root!(, vesicle, tspan, residual; res_options, kwargs...)

Alter the parameters using root-finding with NLsolve.nlsolve. Returns the object returned by NLsolve.nlsolve.

source
VesicleForms.ODEsolveMethod
ODEsolve(vesicle, tspan; plotkws, kwargs...)

Integrate the shape equations from tspan[1] to tspan[end] for the parameters of vesicle.

source
VesicleForms.changeFunction

Changes the value of the variable of a vesicle given as field in steps of dv until target is reached.

change(vesicle, target, tspan)
change(vesicle, target, tspan, field)
change(vesicle, target, tspan, field, dv; silent, stop, store, residual, factor, kwargs...)

Positional arguments

  • vesicle: the vesicle to start with
  • target : the target value to reach
  • tspan : the interval of the arclength to integrate over, typically (0.0, 1.0)
  • field : valid symbol of a property of the vesicle to change (like :v (default) or :m̄)
  • dv : stepsize (default: 1/100 of the distance to target)

Keyword arguments

  • silent : whether to print the progress to stdout (false)
  • stop : callback that can be applied to stop before reaching the target (nothing)
  • store : container to store intermediate results that will be push!ed to it (nothing)
  • residual : function to compute the residual for the non-linear solver (root_residualAV)
  • factor : initial factor for controlling the size of the trust-region
  • additional keyword arguments are passed to ODEshoot_root

Returns the final Vesicle.

source
VesicleForms.matchShotFunction
matchShot(vesicle, tspan)
matchShot(vesicle, tspan, save; dt, alg)

Construct a shape by integrating the shapeEquations in the intervalls [ tspan[1], tspan[end÷2] ] and [ tspan[end], tspan[end÷2] ] with the parameters given by collect( vesicle ). Return ( u, t, prob ).

source
VesicleForms.relax_to_vFunction
relax_to_v(vesicle, target_v, tspan)
relax_to_v(vesicle, target_v, tspan, dv; silent, residual, factor, kwargs...)

Adapt stepwise the reduced volume and adjust parameters via ODEshoot_root! until target_v is reached. Return the final Vesicle.

source
VesicleForms.root_residualAV!Method
root_residualAV!(res, parameters, vesicle, tspan; res_options, kwargs...)

Integrate the shapeEquations using solve_shapeEquations save the difference to endConditionsAV in res.

source
VesicleForms.solve_shapeEquationsMethod

Solve the shape equations with fixed parameters.

solve_shapeEquations(vesicle, tspan; dt, rcut, kwargs...)

Positional arguments

  • vesicle : vesicle to start with
  • tspan : the interval of the arclength to integrate over, typically (0.0, 1.0)

Keyword arguments

  • dt : stepsize for the integration (1e-3)
  • rcut : value of the radial coordinate until which taylor expansions are used instead of integrating the shapeEquations (1e-2)
  • additional keywords are passed to integrate_shapeEquations

Returns (sol::ODESolution, Δu), where sol is the augmented solution object returned from DifferentialEquations.jl and Δu is the difference of the endpoint of the integration and the taylor expansion around the true endpoint evaluated at the actual endpoint.

source

Vesicle prototypes

VesicleForms.prolate_ideal_neck: m̄ = 0.0, v = 0.7071

VesicleForms.noPressureSphere: m̄ = 0.0, v = 1.0

VesicleForms.prolate0_7: m̄ = 0.0, v = 0.7

VesicleForms.prolate0_685: m̄ = 0.0, v = 0.685

VesicleForms.prolate_sqrt2_m_max: m̄ = 1.067, v = 0.7071

VesicleForms.prolate0_9: m̄ = 0.0, v = 0.9

VesicleForms.central_pear: m̄ = 1.152, v = 0.7113

VesicleForms.left_pear: m̄ = 1.22, v = 0.7318

VesicleForms.oblate0_9: m̄ = 0.0, v = 0.9

VesicleForms.prolate0_685m_max: m̄ = 1.199, v = 0.685

VesicleForms.right_pear: m̄ = 1.22, v = 0.7318

VesicleForms.critSphere: m̄ = 0.0, v = 1.0

VesicleForms.oblate0_835: m̄ = 0.0, v = 0.835

VesicleForms.oblate0_65: m̄ = 0.0, v = 0.65

VesicleForms.prolate_v_0_76_m_0_7: m̄ = 0.7, v = 0.76

VesicleForms.oblate0_6: m̄ = 0.0, v = 0.6

VesicleForms.prolate_ideal_neck_m_max: m̄ = 1.206, v = 0.7071

VesicleForms.mexicanHat: m̄ = 0.0, v = 0.4277

Citing VesicleForms.jl

If you use VesicleForms.jl for academic research and wish to cite it, please use the following paper.

Simon Christ, Thomas Litschel, Petra Schwille and Reinhard Lipowsky, Active shape oscillations of giant vesicles with cyclic closure and opening of membrane necks, Soft Matter, 2021, 17, 319

DOI: 10.1039/d0sm00790k

Index