C++ typenameAbstractMesh::AddAppliedForceContribution方法代码示例

void VertexCryptBoundaryForce<DIM>::AddForceContribution(AbstractCellPopulation<DIM>& rCellPopulation)
{
    // Helper variable that is a static cast of the cell population
    VertexBasedCellPopulation<DIM>* p_cell_population = static_cast<VertexBasedCellPopulation<DIM>*>(&rCellPopulation);

    // Throw an exception message if not using a VertexBasedCellPopulation
    if (dynamic_cast<VertexBasedCellPopulation<DIM>*>(&rCellPopulation) == nullptr)
    {
        EXCEPTION("VertexCryptBoundaryForce is to be used with VertexBasedCellPopulations only");
    }

    // Iterate over nodes
    for (typename AbstractMesh<DIM,DIM>::NodeIterator node_iter = p_cell_population->rGetMesh().GetNodeIteratorBegin();
         node_iter != p_cell_population->rGetMesh().GetNodeIteratorEnd();
         ++node_iter)
    {
        double y = node_iter->rGetLocation()[1]; // y-coordinate of node

        // If the node lies below the line y=0, then add the boundary force contribution to the node forces
        if (y < 0.0)
        {
            c_vector<double, DIM> boundary_force = zero_vector<double>(DIM);
            boundary_force[1] = mForceStrength*SmallPow(y, 2);

            node_iter->AddAppliedForceContribution(boundary_force);
        }
    }
}

void TCellDiffusionForce<DIM>::AddForceContribution(AbstractCellPopulation<DIM>& rCellPopulation)
{
    double dt = SimulationTime::Instance()->GetTimeStep();

    // Iterate over the nodes
    for (typename AbstractMesh<DIM, DIM>::NodeIterator node_iter = rCellPopulation.rGetMesh().GetNodeIteratorBegin();
         node_iter != rCellPopulation.rGetMesh().GetNodeIteratorEnd();
         ++node_iter)
    {
        // Get the radius of this node
        unsigned node_index = node_iter->GetIndex();
        double node_radius = node_iter->GetRadius();
        
        // Get cell associated with this index
        CellPtr p_cell = rCellPopulation.GetCellUsingLocationIndex(node_index);

        // Reject if no radius has been set
        if (node_radius == 0.0)
        {
            EXCEPTION("SetRadius() must be called on each Node before calling TCellDiffusionForce::AddForceContribution() to avoid a division by zero error");
        }
        
        //   If the selected cell is a Unlabelled Differentiated T Cell, apply diffusion force contribution.
        if (  (p_cell->GetMutationState()->IsType<TCellMutationState>()) && (p_cell->GetCellProliferativeType()->IsType<DifferentiatedCellProliferativeType>()) 
            && !(p_cell->HasCellProperty<CellLabel>())  )
        {
            double nu = dynamic_cast<AbstractOffLatticeCellPopulation<DIM>*>(&rCellPopulation)->GetDampingConstant(node_index);

            /* Compute the diffusion coefficient D as D = k*T/(6*pi*eta*r), where
             *
             * k = Boltzmann's constant,
             * T = absolute temperature,
             * eta = dynamic viscosity,
             * r = cell radius. */
            double diffusion_const_scaling = GetDiffusionScalingConstant();
            double diffusion_constant = diffusion_const_scaling/node_radius;

            c_vector<double, DIM> force_contribution;
            for (unsigned i=0; i<DIM; i++)
            {
                /* The force on this cell is scaled with the timestep such that when it is
                 * used in the discretised equation of motion for the cell, we obtain the
                 * correct formula
                 *
                 * x_new = x_old + sqrt(2*D*dt)*W
                 *
                 * where W is a standard normal random variable. */
                 
                double xi = RandomNumberGenerator::Instance()->StandardNormalRandomDeviate();
                force_contribution[i] = mStrengthParameter * ((nu*sqrt(2.0*diffusion_constant*dt)/dt)*xi);
            }
            node_iter->AddAppliedForceContribution(force_contribution);
        }
    }
    
    
}

void MeshBasedCellPopulationWithGhostNodes<DIM>::ApplyGhostForces(){

    // Initialise vector of forces on ghost nodes
    std::vector<c_vector<double, DIM> > drdt(this->GetNumNodes());
    for (unsigned i=0; i<drdt.size(); i++)
    {
        drdt[i] = zero_vector<double>(DIM);
    }

    // Calculate forces on ghost nodes
    for (typename MutableMesh<DIM, DIM>::EdgeIterator edge_iterator = static_cast<MutableMesh<DIM, DIM>&>((this->mrMesh)).EdgesBegin();
        edge_iterator != static_cast<MutableMesh<DIM, DIM>&>((this->mrMesh)).EdgesEnd();
        ++edge_iterator)
    {
        unsigned nodeA_global_index = edge_iterator.GetNodeA()->GetIndex();
        unsigned nodeB_global_index = edge_iterator.GetNodeB()->GetIndex();

        c_vector<double, DIM> force = CalculateForceBetweenGhostNodes(nodeA_global_index, nodeB_global_index);

        if (!this->mIsGhostNode[nodeA_global_index])
        {
            drdt[nodeB_global_index] -= force;
        }
        else
        {
            drdt[nodeA_global_index] += force;

            if (this->mIsGhostNode[nodeB_global_index])
            {
                drdt[nodeB_global_index] -= force;
            }
        }
    }

    for (typename AbstractMesh<DIM,DIM>::NodeIterator node_iter = this->mrMesh.GetNodeIteratorBegin();
         node_iter != this->mrMesh.GetNodeIteratorEnd();
         ++node_iter)
    {
        unsigned node_index = node_iter->GetIndex();
        if (this->mIsGhostNode[node_index])
        {
            node_iter->ClearAppliedForce();
            node_iter->AddAppliedForceContribution(drdt[node_index]);
        }
    }

};