hamiltonian/Maxwell.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 "Maxwell.h"

#include "HHelper.h"
#include "bio/Cell.h"
#include "sim/CoreData.h"

#include <cassert>

using namespace std;
using namespace boost::property_tree;
using namespace SimPT_Sim;

namespace SimPT_Default {
namespace Hamiltonian {

Maxwell::Maxwell(const CoreData& cd)
{
        assert( cd.Check() && "CoreData not ok in Maxwell Hamiltonian");
        Initialize(cd);
}

void Maxwell::Initialize(const CoreData& cd)
{
        const auto& p = cd.m_parameters->get_child("cell_mechanics");

        m_elastic_modulus                = p.get<double>("elastic_modulus");
        m_lambda_bend                    = p.get<double>("lambda_bend");
        m_lambda_cell_length             = p.get<double>("lambda_celllength");
        m_lambda_length                  = p.get<double>("lambda_length");
        m_rp_stiffness                   = p.get<double>("relative_perimeter_stiffness");
        m_target_node_distance           = p.get<double>("target_node_distance");
}

double Maxwell::operator()(Cell* cell)
{
        double h = 0.0;

        if (!cell->IsBoundaryPolygon()) {

                // --------------------------------------------------------------------------------------------
                // Cell pressure:
                // --------------------------------------------------------------------------------------------
                h += -cell->GetSolute() * log(cell->GetArea());

                // --------------------------------------------------------------------------------------------
                // Elastic wall:
                // --------------------------------------------------------------------------------------------
                h += m_elastic_modulus * HHelper::ElasticWallTerm(cell, m_rp_stiffness);

                // --------------------------------------------------------------------------------------------
                // Edge length:
                // --------------------------------------------------------------------------------------------
                if (abs(m_lambda_length) > 1.0e-7) {
                        h += m_lambda_length
                                * HHelper::EdgeTerm(cell, m_rp_stiffness, m_target_node_distance);
                }

                // --------------------------------------------------------------------------------------------
                // Cell length constraint:
                // --------------------------------------------------------------------------------------------
                if (abs(m_lambda_cell_length) > 1.0e-7) {
                        h += m_lambda_cell_length * HHelper::CellLengthTerm(cell);
                }
        }

        // --------------------------------------------------------------------------------------------
        // Vertex bending constraint
        // --------------------------------------------------------------------------------------------
        if (abs(m_lambda_bend) > 1.0e-7) {
                h +=  m_lambda_bend * HHelper::BendingTerm(cell);
        }

        return h;
}

} // namespace
} // namespace