The Source File for the Action Class Virus
¶
At the beginning of Virus.cpp
we must put the header files we need:
#include <omp.h>
#include "MessLoggerT.h"
#include "clsutils.cpp"
#include "ParamProvider2.h"
#include "SPopulation.h"
#include "SCell.h"
#include "WELL512.h"
#include "QDFUtils.h"
#include "QDFUtilsT.h"
#include "Virus.h"
We also need to define attribute names to read and write the action’s attributes.
template<typename T>
const char *Virus<T>::asNames[] = {
ATTR_VIRUS_INFECTION_PROB_NAME,
ATTR_VIRUS_INITIAL_LOAD_NAME,
ATTR_VIRUS_GROWTH_RATE_NAME,
ATTR_VIRUS_CONTAGION_LEVEL_NAME,
ATTR_VIRUS_LETHALITY_LEVEL_NAME,
};
constructor
¶
There is not much to do in the constructor: intialisation of the variables, and saving the attribute names.
template<typename T>
Virus<T>::Virus(SPopulation<T> *pPop, SCellGrid *pCG, std::string sID, WELL512 **apWELL)
: Action<T>(pPop,pCG,ATTR_VIRUS_NAME,sID),
m_apWELL(apWELL),
m_fInfectionProb(0),
m_fInitialLoad(0),
m_fGrowthRate(0),
m_fContagionLevel(0),
m_fLethalityLevel(0),
m_iNumThreads(omp_get_max_threads()) {
this->m_vNames.insert(this->m_vNames.end(), asNames, asNames+sizeof(asNames)/sizeof(char*));
// create array with one int map for each thread
m_amCellAgents = new cellagentmap[m_iNumThreads];
// create array with one map for each thread
m_amAgentInfects = new agentloadmap[m_iNumThreads];
}
destructor
¶
The destructor deletes the arrays that were allocated in the constructor
template<typename T>
Virus<T>::~Virus() {
// delete the arrays
delete[] m_amCellAgents;
delete[] m_amAgentInfects;
}
initialize
¶
The method initialize is called before execute()
.
Here the array of maps relating cell indexes to agent indexes (m_amCellAGents
) is filled,
and the array of maps relating agent IDs to virus loads (m_amAgentInfects
) is resaet.
Note that since initialize
is not called in a parllelized context,
we can use parallelized loops with the full set of available threads.
template<typename T>
int Virus<T>::initialize(float fT) {
int iResult = 0;
// clear the agent arrays and the infection registerr
for (int iT = 0; iT < m_iNumThreads; iT++) {
m_amCellAgents[iT].clear();
m_amAgentInfects[iT].clear();
}
int iFirstAgent = this->m_pPop->getFirstAgentIndex();
int iLastAgent = this->m_pPop->getLastAgentIndex();
#pragma omp parallel for
for (int iA = iFirstAgent; iA <= iLastAgent; iA++) {
int iT = omp_get_thread_num();
T* pA = &(this->m_pPop->m_aAgents[iA]);
m_amCellAgents[iT][pA->m_iCellIndex]->push_back(iA);
}
return iResult;
}
execute
¶
The execute()
method is applied to every agent (in parallel).
This method is usually the one most of the action’s work is done.
Here is what we need the virus to do for each agent:
- let the the agent’s virus load grow logistically
- if the agent’s virus load exceeds the lethality level, register it for death
- if the agent’s virus load exceeds the contagion level, it infects all agents in the cell
template<typename T>
int Virus<T>::execute(int iAgentIndex, float fT) {
int iResult = 0;
int iT = omp_get_thread_num();
T *pa = &(this->m_pPop->m_aAgents[iAgentIndex]);
// we only handle agents not already marked as dead
if (pa->m_iLifeState > 0) {
// perform the logistic growth of the agent's viral load
float fM = pa->m_fViralLoad;
pa->m_fViralLoad += m_fGrowthRate*fM*(1-fM);
int iCellIndex = pa->m_iCellIndex;
// check if the agent must die
if (pa->m_fViralLoad > m_fLethalityLevel) {
this->m_pPop->registerDeath(iCellIndex, iAgentIndex);
// is the agent contagious?
} else if (pa->m_fViralLoad > m_fContagionLevel) {
// loop through agents in cell
for (uint i = 0; i < m_amCellAgents[iT][iCellIndex].size(); i++) {
int iOtherAgent = m_amCellAgents[iCellIndex][i];
// we only infect other agents
if (iOtherAgent != iAgentIndex) {
T *pao = &(this->m_pPop->m_aAgents[iOtherAgent]);
// infect!
double dR = this->m_apWELL[omp_get_thread_num()]->wrandd();
// we reduce the effective infection probability by multiplying it with (1-immunity)
if (dR < m_fInfectionProb * (1-pao->m_fImmunity)) {
// register infection load for agent
m_amAgentInfects[iT][iOtherAgent] += m_fInitialLoad;
}
}
}
}
}
return iResult;
}
finalize
¶
The method finalize
is called after all calls to execute()`.
Here we finally add to each agent’s virus load the colllected virus load
passed from its neighbors. During execute()
several maps may have entries
for the same agent.
template<typename T>
int Virus<T>::finalize(float fT) {
int iResult = 0;
for (int iT = 0; iT < m_iNumThreads; iT++) {
for (agentloadmap::const_iterator it = m_amAgentInfects[iT].begin(); it != m_amAgentInfects[iT].end(); ++it) {
// first: agent index
// second: infecttload
T *pa = &(this->m_pPop->m_aAgents[it->first]);
pa->m_fViralLoad += it->second;
}
}
return iResult;
}
The remaining methods deal with reading and writing the virus’ attributes.
extractAttributesQDF
¶
This method reads virus attributes from a QDF file.
template<typename T>
int Virus<T>::extractAttributesQDF(hid_t hSpeciesGroup) {
int iResult = 0;
if (iResult == 0) {
iResult = qdf_extractAttribute(hSpeciesGroup, ATTR_VIRUS_INFECTION_PROB_NAME, 1, &m_fInfectionProb);
if (iResult != 0) {
LOG_ERROR("[Virus] couldn't read attribute [%s]", ATTR_VIRUS_INFECTION_PROB_NAME);
}
}
if (iResult == 0) {
iResult = qdf_extractAttribute(hSpeciesGroup, ATTR_VIRUS_INITIAL_LOAD_NAME, 1, &m_fInitialLoad);
if (iResult != 0) {
LOG_ERROR("[Virus] couldn't read attribute [%s]", ATTR_VIRUS_INITIAL_LOAD_NAME);
}
}
if (iResult == 0) {
iResult = qdf_extractAttribute(hSpeciesGroup, ATTR_VIRUS_GROWTH_RATE_NAME, 1, &m_fGrowthRate);
if (iResult != 0) {
LOG_ERROR("[Virus] couldn't read attribute [%s]", ATTR_VIRUS_GROWTH_RATE_NAME);
}
}
if (iResult == 0) {
iResult = qdf_extractAttribute(hSpeciesGroup, ATTR_VIRUS_CONTAGION_LEVEL_NAME, 1, &m_fContagionLevel);
if (iResult != 0) {
LOG_ERROR("[Virus] couldn't read attribute [%s]", ATTR_VIRUS_CONTAGION_LEVEL_NAME);
}
}
if (iResult == 0) {
iResult = qdf_extractAttribute(hSpeciesGroup, ATTR_VIRUS_LETHALITY_LEVEL_NAME, 1, &m_fLethalityLevel);
if (iResult != 0) {
LOG_ERROR("[Virus] couldn't read attribute [%s]", ATTR_VIRUS_LETHALITY_LEVEL_NAME);
}
}
return iResult;
}
writeAttributesQDF
¶
This method writes the virus’ attributes to a QDF file.
template<typename T>
int Virus<T>:: writeAttributesQDF(hid_t hSpeciesGroup) {
int iResult = 0;
iResult += qdf_insertAttribute(hSpeciesGroup, ATTR_VIRUS_INFECTION_PROB_NAME, 1, &m_fInfectionProb);
iResult += qdf_insertAttribute(hSpeciesGroup, ATTR_VIRUS_INITIAL_LOAD_NAME, 1, &m_fInitialLoad);
iResult += qdf_insertAttribute(hSpeciesGroup, ATTR_VIRUS_GROWTH_RATE_NAME, 1, &m_fGrowthRate);
iResult += qdf_insertAttribute(hSpeciesGroup, ATTR_VIRUS_CONTAGION_LEVEL_NAME, 1, &m_fContagionLevel);
iResult += qdf_insertAttribute(hSpeciesGroup, ATTR_VIRUS_LETHALITY_LEVEL_NAME, 1, &m_fLethalityLevel);
return iResult;
}
tryGetAttributes
¶
This method reads the virus attributes passed via the xml population file<pop_xml_ref>.
template<typename T>
int Virus<T>::tryGetAttributes(const ModuleComplex *pMC) {
int iResult = -1;
const stringmap &mParams = pMC->getAttributes();
if (this->checkAttributes(mParams) == 0) {
iResult = 0;
iResult += getAttributeVal(mParams, ATTR_VIRUS_INFECTION_PROB_NAME, &m_fInfectionProb);
iResult += getAttributeVal(mParams, ATTR_VIRUS_INITIAL_LOAD_NAME, &m_fInitialLoad);
iResult += getAttributeVal(mParams, ATTR_VIRUS_GROWTH_RATE_NAME, &m_fGrowthRate);
iResult += getAttributeVal(mParams, ATTR_VIRUS_CONTAGION_LEVEL_NAME, &m_fContagionLevel);
iResult += getAttributeVal(mParams, ATTR_VIRUS_LETHALITY_LEVEL_NAME, &m_fLethalityLevel);
}
return iResult;
}
The entire file can be found here: Virus.cpp
Implementation of the population class VirusHostPop
¶
Since our Virus` action makes certain assumptions about the agents, we can’t use an arbitary population class. In particular, it needs fields for the virus load and immunity.
The header file¶
First we have the includes for the header files of all required
actions (among them Virus.h
), and the header file for the
population’s base class SPopulation
#include "GetOld.h"
#include "ATanDeath.h"
#include "RandomMove.h"
#include "Fertility.h"
#include "Verhulst.h"
#include "RandomPair.h"
#include "Virus.h"
#include "SPopulation.h"
Some constants needed to extract mutation rate and inheritance types
const std::string VAR_VIRUSHOST_MUT_RATE_NAME = "MutationRate";
const std::string VAR_VIRUSHOST_IMM_INHERIT_NAME = "ImmunityInheritance";
const std::string INH_TYPES[] = {
"mix",
"mat",
"pat",
"min",
"max",
};
const int INH_MIX = 0;
const int INH_MAT = 1;
const int INH_PAT = 2;
const int INH_MIN = 3;
const int INH_MAX = 4;
Now the declaration of the agent structure: we have members for agent fertility, and the members for the viral infections virus load, immunity and incoming virus load.
struct VirusHostAgent : Agent {
float m_fAge;
float m_fLastBirth;
int m_iMateIndex;
float m_fViralLoad;
float m_fImmunity;
};
The declaration of the VirusHostPop
class derived from SPopulation
:
class VirusHostPop : public SPopulation<VirusHostAgent> {
The public methods are the usual: constructor, destructor and methods to add agent data from a populatiom dat file<pop_dat_ref> and to write agent data to a QDF file. In addition to that there is a method for setting the inital state of the baby (which we’ve also encountered before), and a method to get a parameter for the population from the populaition xml file<pop_xml_ref> (in this case we need the mutation rate, and the inheritance type).
public:
VirusHostPop(SCellGrid *pCG, PopFinder *pPopFinder, int iLayerSize, IDGen **apIDG, uint32_t *aulState, uint *aiSeeds);
virtual ~VirusHostPop_SexualPop();
int addPopSpecificAgentData(int iAgentIndex, char **ppData);
void addPopSpecificAgentDataTypeQDF(hid_t *hAgentDataType);
int makePopSpecificOffspring(int iAgent, int iMother, int iFather);
int getPopParams(const stringmap &mVarDefs);
In the protected
section there are pointers to the actions of this population
and the variable for the mutation rate.
protected:
GetOld<VirusHostAgent> *m_pGO;
ATanDeath<VirusHostAgent> *m_pAD;
RandomMove<VirusHostAgent> *m_pRM;
Fertility<VirusHostAgen> *m_pFert;
Verhulst<VirusHostAgen> *m_pVerhulst;
RandomPair<VirusHostAgen> *m_pPair;
Virus<VirusHostAgen> *m_pVirus;
float m_fMutationRate;
constructor¶
Here we instantiate all actions and initialize the fields.
VirusHostPop::VirusHostPop(SCellGrid *pCG, PopFinder *pPopFinder, int iLayerSize, IDGen **apIDG, uint32_t *aulState, uint *aiSeeds)
: SPopulation<VirusHostAgent>(pCG, pPopFinder, iLayerSize, apIDG, aulState, aiSeeds) {
m_pGO = new GetOld<VirusHostAgent>(this, m_pCG, "");
m_pAD = new ATanDeath<VirusHostAgent>(this, m_pCG, "", m_apWELL);
m_pRM = new RandomMove<VirusHostAgent>(this, m_pCG, "", m_apWELL);
m_pFert = new Fertility<VirusHostAgent>(this, m_pCG, "");
m_pVerhulst = new Verhulst<VirusHostAgent>(this, m_pCG, "", m_apWELL);
m_pPair = new RandomPair<VirusHostAgent>(this, m_pCG, "", m_apWELL);
m_pVirus = new Virus<VirusHostAgent>(this, m_pCG, "", m_apWELL);
m_prio.addAction(m_pGO);
m_prio.addAction(m_pAD);
m_prio.addAction(m_pRM);
m_prio.addAction(m_pFert);
m_prio.addAction(m_pVerhulst);
m_prio.addAction(m_pPair);
m_prio.addAction(m_pVirus);
m_fMutationRate = 0.0;
m_sInheritType = "mix";
}
destructor¶
The destructor deletes all action created in the constructor.
VirusHostPop::~VirusHostPop() {
if (m_pGO != NULL) {
delete m_pGO;
}
if (m_pAD != NULL) {
delete m_pAD;
}
if (m_pRM != NULL) {
delete m_pRM;
}
if (m_pFert != NULL) {
delete m_pFert;
}
if (m_pVerhulst != NULL) {
delete m_pVerhulst;
}
if (m_pPair != NULL) {
delete m_pPair;
}
if (m_pVirus != NULL) {
delete m_pVirus;
}
}
getPopParams¶
This method tries to extract mutation rate and inheritance type from the stringmap created by the population XML reader.
int VirusHostPop::getPopParams(const stringmap &mVarDefs) {
int iResult = -1;
stringmap::const_iterator it = mVarDefs.find(VAR_VIRUSHOST_MUT_RATE_NAME);
if (it != mVarDefs.end()) {
float fMutTemp = 0;
if (strToNum(it->second, &fMutTemp)) {
m_fMutationRate = fMutTemp;
printf("[VirusHostPop::getPopParams] have mutation rate %f\n", m_fMutationRate);
iResult = 0;
} else {
// couldn't convert string to float
}
} else {
// there is no key VAR_VIRUSHOST_MUT_RATE_NAME
}
it = mVarDefs.find(VAR_VIRUSHOST_IMM_INHERIT_NAME);
if (it != mVarDefs.end()) {
std::string sInheritType = it->second;
// handle the "+"
if (!sInheritType.ends_with("+")) {
m_fMutationRate = 0;
} else {
// drop the plus
sInheritType = sInheritType.substr(0, sInheritType.size()-1);
}
// check if valid inherit type
bool bSearching = true;
for (int i = 0; bSearching && (i < 5); i++) {
if (sInheritType == INH_TYPES[i]) {
bSearching = false;
}
}
if (bSearching) {
printf("[VirusHostPop::getPopParams] unknown inherit type [%s]\n", sInheritType.c_str());
} else {
m_sInheritType = sInheritType;
printf("[VirusHostPop::getPopParams] have inherit type [%s]\n", m_sInheritType.c_str());
iResult = 0;
}
} else {
// there is no key VAR_VIRUSHOST_IMM_INHERIT_NAME
}
return iResult;
}
makePopSpecificOffspring¶
This method is called when a baby is born in order to initialize it’s field. In particular, the baby’s immunity level is determined by the inhertiance type.
int VirusHostPop::makePopSpecificOffspring(int iAgent, int iMother, int iFather) {
int iResult = 0;
m_aAgents[iAgent].m_fAge = 0.0;
m_aAgents[iAgent].m_fLastBirth = 0.0;
m_aAgents[iAgent].m_iMateIndex = -3;
// SPopulation assigns random genders to a baby, ehich is ok here
m_aAgents[iAgent].m_fViralLoad = 0;
// register this date as the last birth of the mother
m_aAgents[iAgent].m_fLastBirth = m_fCurTime;
// the offspring's immunity is the average of the parents' immunities plus a little fuzz factor.
float fImm = 0.0;
if (m_sInheritType == INH_TYPES[INH_MIX]) {
fImm = (m_aAgents[iMother].m_fImmunity + m_aAgents[iFather].m_fImmunity)/2;
} else if (m_sInheritType == INH_TYPES[INH_MAT]) {
fImm = m_aAgents[iMother].m_fImmunity;
} else if (m_sInheritType == INH_TYPES[INH_PAT]) {
fImm = m_aAgents[iFather].m_fImmunity;
} else if (m_sInheritType == INH_TYPES[INH_MIN]) {
fImm = (m_aAgents[iMother].m_fImmunity < m_aAgents[iFather].m_fImmunity)?m_aAgents[iMother].m_fImmunity:m_aAgents[iFather].m_fImmunity;
} else if (m_sInheritType == INH_TYPES[INH_MAX]) {
fImm = (m_aAgents[iMother].m_fImmunity > m_aAgents[iFather].m_fImmunity)?m_aAgents[iMother].m_fImmunity:m_aAgents[iFather].m_fImmunity;
} else {
// shouldn't happen
iResult = -1;
}
//mutate
fImm += m_apWELL[omp_get_thread_num()]->wgauss(m_fMutationRate);
if (fImm < 0) {
fImm = 0;
}else if (fImm > 1) {
fImm = 1;
}
m_aAgents[iAgent].m_fImmunity = fImm;
if (m_aAgents[iMother].m_fViralLoad > 0) {
m_aAgents[iAgent].m_fViralLoad = 0.2;
printf("Babyinfect: %d -> %d\nn", iMother, iAgent);
}
//printf("[VirusHostPop::makePopSpecificOffsprin] baby %d vl %f, imm %f\n", iAgent, m_aAgents[iAgent].m_fViralLoad, m_aAgents[iAgent].m_fImmunity);
return iResult;
}
addPopSpecificAgentData¶
This method reads agent specific data from the population dat file<pop_dat_ref>.
int VirusHostPop::addPopSpecificAgentData(int iAgentIndex, char **ppData) {
int iResult = 0;
if (iResult == 0) {
iResult = addAgentDataSingle(ppData, &m_aAgents[iAgentIndex].m_fAge);
}
if (iResult == 0) {
iResult = addAgentDataSingle(ppData, &m_aAgents[iAgentIndex].m_fLastBirth);
}
if (iResult == 0) {
iResult = addAgentDataSingle(ppData, &m_aAgents[iAgentIndex].m_fViralLoad);
}
if (iResult == 0) {
iResult = addAgentDataSingle(ppData, &m_aAgents[iAgentIndex].m_fImmunity);
}
return iResult;
}
addPopSpecificAgentDataTypeQDF¶
This method creates a HDF data type to handlle I/O of agents.
void VirusHostPop::addPopSpecificAgentDataTypeQDF(hid_t *hAgentDataType) {
VirusHostAgent va;
H5Tinsert(*hAgentDataType, "Age", qoffsetof(va, m_fAge), H5T_NATIVE_FLOAT);
H5Tinsert(*hAgentDataType, "LastBirth", qoffsetof(va, m_fLastBirth), H5T_NATIVE_FLOAT);
H5Tinsert(*hAgentDataType, "ViralLoad", qoffsetof(va, m_fViralLoad), H5T_NATIVE_FLOAT);
H5Tinsert(*hAgentDataType, "Immunity", qoffsetof(va, m_fImmunity), H5T_NATIVE_FLOAT);
}