cell_housekeep/SmithPhyllotaxis.cpp

Go to the documentation of this file.

/*
 * Copyright 2011-2016 Universiteit Antwerpen
 *
 * Licensed under the EUPL, Version 1.1 or  as soon they will be approved by
 * the European Commission - subsequent versions of the EUPL (the "Licence");
 * You may not use this work except in compliance with the Licence.
 * You may obtain a copy of the Licence at: http://ec.europa.eu/idabc/eupl5
 *
 * Unless required by applicable law or agreed to in writing, software
 * distributed under the Licence is distributed on an "AS IS" basis,
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 * See the Licence for the specific language governing
 * permissions and limitations under the Licence.
 */
#include "SmithPhyllotaxis.h"


#include "bio/Cell.h"
#include "bio/Mesh.h"

#include <cmath>

namespace SimPT_Default {
namespace CellHousekeep {

using namespace std;
using namespace boost::property_tree;

SmithPhyllotaxis::SmithPhyllotaxis(const CoreData& cd)
{
        Initialize(cd);
}

void SmithPhyllotaxis::Initialize(const CoreData& cd)
{
        m_cd = cd;
        auto& p = m_cd.m_parameters;

        // Initialize model parameters from ptree
        m_cell_expansion_rate = p->get<double>("cell_mechanics.cell_expansion_rate");
        m_elastic_modulus = p->get<double>("cell_mechanics.elastic_modulus");
        m_response_time = p->get<double>("cell_mechanics.response_time");
        m_time_step = p->get<double>("model.time_step");
        m_viscosity_const = p->get<double>("cell_mechanics.viscosity_const");
        // Calculate additional parameters
        m_area_incr = m_cell_expansion_rate * m_time_step;
}

void SmithPhyllotaxis::operator()(Cell* cell)
{
        const string ham_select = m_cd.m_parameters->get<string> ("model.mc_hamiltonian");

        const double IAA = cell->GetChemical(0);
        const double t_area = cell->GetTargetArea();
        const double t_length = cell->GetTargetLength();

        // Logistic function centered around IAA concentration 5 seems to give more
        // localized growth than the expression used in AuxinGrowth
        //const double update_t_area    = t_area + 1.0 / (1.0 + exp(-IAA + 5.0)) * m_area_incr;
        // TODO: reverted for now, DDV suggested it's more adopted in the field.
        const double update_t_area = t_area + IAA / (1.0 + IAA) * m_area_incr;
        // The update in cell length (longest axis of the cell ellipse) is chosen
        // such that l / sqrt(A) is a constant which means that the cell should be
        // as "square" (or round) as possible.
        const double update_t_length = t_length * sqrt(update_t_area / t_area);

        cell->SetTargetArea(update_t_area);
        cell->SetTargetLength(update_t_length);

        if (ham_select == "ElasticWall") {
                // Update the rest length of wall
                // TODO: 1.50 (maximal extension) and 1.25 (irreversible extension) are hidden parameters
                for (list<Wall*>::iterator i = cell->GetWalls().begin(); i != cell->GetWalls().end(); ++i) {
                        Wall* w = *i;
                        if (w->GetLength() > 1.50 * w->GetRestLength()) {
                                // Irreversible extension
                                w->SetRestLength(w->GetLength() / 1.25);
                        }
                }
        }

        if (ham_select == "Maxwell") {
                // Update the solute of cell and the rest length of wall
                const double solute = m_viscosity_const * sqrt(cell->GetArea());
                const double update_solute = solute
                        - (solute - cell->GetSolute()) * exp(-0.1 * m_time_step / m_response_time);
                cell->SetSolute(update_solute);

                for (list<Wall*>::iterator i = cell->GetWalls().begin(); i != cell->GetWalls().end(); ++i) {
                        Wall* w = *i;

                        const double update_rest_length = w->GetLength()
                                - (w->GetLength() - w->GetRestLength())
                                        * exp(-0.01 * m_elastic_modulus * m_time_step / m_viscosity_const);
                        w->SetRestLength(update_rest_length);
                }
        }
}

} // namespace
} // namespace