Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F82303783
polarSignaling.cpp
No One
Temporary
Actions
Download File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Subscribers
None
File Metadata
Details
File Info
Storage
Attached
Created
Tue, Sep 10, 16:15
Size
4 KB
Mime Type
text/x-c
Expires
Thu, Sep 12, 16:15 (1 d, 21 h)
Engine
blob
Format
Raw Data
Handle
20676819
Attached To
R9411 tisue modeling
polarSignaling.cpp
View Options
#include "polarSignaling.hpp"
using namespace std;
PolarSignaling::PolarSignaling(){
tau_pol = 0.01;
}
void PolarSignaling::operator()(const state_type &y, state_type &dydt,double t){
// only allow exchange with neighboring boundaries
// position-dependent pull factor
// cout<<"starting us"<<endl;
updateSignal(y);
// cout<<"uss done"<<endl;
int n= lattice->getNbCells();
for (int c=0;c<n;c++){
Cell& ce = lattice->getCell(c);
if(!ce.isEmpty() && !ce.isFrozen()){
double sum=0.0;
BoundaryIterator it(&ce);
for(Boundary *b=it.first();!it.isLast();b=it.next()){
// cout<<"idx 0 "<<b->id()<<endl;
uint idx = b->stateVectorIndex();
uint cidx = b->cell()->id();
if(b->oppositeCell()){
uint oidx = b->oppositeCell()->id();
// cout<<"idx "<<idx<<" "<<oidx<<endl;
// cout<<"idx2 "<<idx<<" "<<oidx<<endl;
Vec2 v = b->positionRelativeToCell();
try{v.Normalize();}
catch(...){v=Vec2(0,0);}
// cout<<"position "<<v<<endl;
// cout<<"vi "<<v<<endl;
// if (isnan(y[idx])) {cout<<"cell "<<c<<" vec "<<v<<endl;exit(0);}
// v *= fabs(y[idx])/v.Norm();
v *= fabs(y[idx]);///v.Norm(); //why not 1/v.Norm()???
double dv = tau_pol*sin(v.Dot(m_cellPolarity[oidx]+m_cellPolarity[cidx]));
// cout<<"dv "<<dv<<" / "<<v.Dot(m_cellPolarity[oidx]) << ":"<< m_cellPolarity[oidx] << " ++ " << v<<endl;
dydt[idx] = dv;
sum += dv;
// cout<<"dv out "<<sum<<endl;
}
}
sum /= ce.getNbBoundaries();
// cout<<"sum "<<ce.id()<<": "<<sum<<endl;
// remove average dydt
for(Boundary *b=it.first();!it.isLast();b=it.next()){
uint idx = b->stateVectorIndex();
dydt[idx] -= sum;
}
}
}
}
/*
void PolarSignaling::operator()(const state_type &y, state_type &dydt,double t){
// only allow exchange with neighboring boundaries
// position-dependent pull factor
int n= lattice->getNbCells();
for (int c=0;c<n;c++){
Cell& ce = lattice->getCell(c);
if(!ce.isFrozen()){
BoundaryIterator it(&ce);
for(Boundary *b=it.first();!it.isLast();b=it.next()){
uint idx = b->stateVectorIndex();
uint oidx = b->opposite()->stateVectorIndex();
uint pidx = b->prev()->stateVectorIndex();
uint nidx = b->next()->stateVectorIndex();
uint opidx = b->prev()->opposite()->stateVectorIndex();
uint onidx = b->next()->opposite()->stateVectorIndex();
Vec2 v = b->positionRelativeToCell();
Vec2 pv = b->prev()->positionRelativeToCell();
Vec2 nv = b->next()->positionRelativeToCell();
double pfac = fabs(fabs(v.sinAngle(Vec2(1,0)))-fabs(pv.sinAngle(Vec2(1,0))));
double nfac = fabs(fabs(v.sinAngle(Vec2(1,0)))-fabs(nv.sinAngle(Vec2(1,0))));
dydt[idx] = tau_pol*y[idx]*
(pfac*y[pidx]*(y[oidx]-y[opidx])+nfac*y[nidx]*(y[oidx]-y[onidx]));
}
}
}
}
*/
float PolarSignaling::getWeight(Boundary *b){
return 1.0;
}
void PolarSignaling::initBoundary(Boundary *bo, state_type *state){
// compute position and share across boundaries
uint idx = bo->stateVectorIndex();
uint n= bo->cell()->getNbBoundaries();
double pol = (RND(0.001)+0.9995)/n;
(*state)[idx]=pol;
m_polarity[bo->id()]=pol;
cout<<pol<<endl;
}
void PolarSignaling::updateSignal(const state_type& state){
uint nbcells = lattice->getNbCells();
for (uint c=0;c<nbcells;c++){
Cell& ce = lattice->getCell(c);
double sum =0.0;
m_cellPolarity[ce.id()]=Vec2(0,0);
if(!ce.isEmpty() && ! ce.isFrozen()){
BoundaryIterator it(&ce);
for(Boundary *b=it.first();!it.isLast();b=it.next()){
double vi = state[b->stateVectorIndex()];
m_polarity[b->id()] =vi;
Vec2 v = b->positionRelativeToCell(true);
// cout<<"us "<< vi<< " : "<<v<<endl;
try{
v.Normalize();
}catch(...){v=Vec2(0,0);}
m_cellPolarity[ce.id()] +=v*vi;
sum += vi;
}
if(sum>EPSILON) m_cellPolarity[ce.id()] /=sum;
}
}
// cout<<"us done"<<endl;
}
void PolarSignaling::initState(state_type& state, double *signal){
uint nbcells = lattice->getNbCells();
for (uint c=0;c<nbcells;c++){
Cell& ce = lattice->getCell(c);
if(!ce.isEmpty()){
BoundaryIterator it(&ce);
for(Boundary *b=it.first();!it.isLast();b=it.next()){
m_polarity[b->id()] = state[b->stateVectorIndex()]=signal[b->id()];
}
}
}
}
void PolarSignaling::initState(state_type& state){
CellInteraction::initState(state);
// uint nbcells = lattice->getNbCells();
// for (uint c=0;c<nbcells;c++){
// Cell& ce = lattice->getCell(c);
// if(!ce.isEmpty()){
// BoundaryIterator it(&ce);
// for(Boundary *b=it.first();!it.isLast();b=it.next()){
// state[b->stateVectorIndex()]=m_polarity[b->id()];
// }
// }
// }
}
void PolarSignaling::fillVector(double *v, int i)const{
uint cnt=0;
switch(i){
case 0:
for(uint i=0;i<m_polarity.size();i++){
v[i]=m_polarity[i];
}
break;
case 1:
for(uint i=0;i<m_cellPolarity.size();i++){
v[cnt++]=m_cellPolarity[i][0];
v[cnt++]=m_cellPolarity[i][1];
}
break;
default:break;
}
}
Event Timeline
Log In to Comment