Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F87479208
collierSignaling.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
Sat, Oct 12, 22:15
Size
4 KB
Mime Type
text/x-c
Expires
Mon, Oct 14, 22:15 (1 d, 21 h)
Engine
blob
Format
Raw Data
Handle
21605536
Attached To
R9411 tisue modeling
collierSignaling.cpp
View Options
#ifndef __COLLIER_SIGNALING_CPP__
#define __COLLIER_SIGNALING_CPP__
#include "collierSignaling.hpp"
#include <assert.h>
using
namespace
std
;
//#define ADIABATIC
template
<
class
Coupling
>
CollierSignaling
<
Coupling
>::
CollierSignaling
(){
v
=
1.0
;
coupling
.
init
();
// set_hill_params(1.0,0.2,3.0);
}
// notch: y[cell()->stateVectorIndex()]
// delta: y[cell()->stateVectorIndex()+1]
template
<
class
Coupling
>
void
CollierSignaling
<
Coupling
>::
operator
()(
const
state_type
&
y
,
state_type
&
dydt
,
double
t
){
updateSignal
(
y
);
int
n
=
lattice
->
getNbCells
();
for
(
int
c
=
0
;
c
<
n
;
c
++
){
Cell
&
ce
=
lattice
->
getCell
(
c
);
if
(
!
ce
.
isEmpty
()
&&
!
ce
.
isFrozen
()){
double
ndelta
=
0.0
;
double
neighcnt
=
0.0
;
double
fac
=
1.0
;
BoundaryIterator
it
(
&
ce
);
for
(
Boundary
*
b
=
it
.
first
();
!
it
.
isLast
();
b
=
it
.
next
()){
if
(
b
->
oppositeCell
()
&&
!
b
->
oppositeCell
()
->
isFrozen
()){
fac
=
b
->
length
()
/
sqrt
(
b
->
oppositeCell
()
->
perimeter
());
#ifndef ADIABATIC
ndelta
+=
MAX
(
0.0
,
y
[
b
->
oppositeCell
()
->
stateVectorIndex
()
+
1
])
*
fac
;
#else
ndelta
+=
m_delta
[
b
->
oppositeCell
()
->
id
()]
*
fac
;
#endif
// neighcnt++;
neighcnt
+=
fac
;
}
}
if
(
neighcnt
>
0
)
ndelta
/=
neighcnt
;
uint
idx
=
ce
.
stateVectorIndex
();
dydt
[
idx
]
=
coupling
.
f_func
(
ndelta
)
-
y
[
idx
];
assert
(
!
isnan
(
ndelta
));
assert
(
!
(
isnan
(
coupling
.
f_func
(
ndelta
))
&&
cout
<<
ndelta
));
assert
(
!
isnan
(
dydt
[
idx
]));
#ifndef ADIABATIC
dydt
[
idx
+
1
]
=
v
*
(
coupling
.
g_func
(
y
[
idx
])
-
y
[
idx
+
1
]);
assert
(
!
isnan
(
dydt
[
idx
+
1
]));
#endif
}
}
}
template
<
class
Coupling
>
void
CollierSignaling
<
Coupling
>::
updateSignal
(
const
state_type
&
y
){
#ifdef ADIABATIC
int
n
=
lattice
->
getNbCells
();
for
(
int
c
=
0
;
c
<
n
;
c
++
){
Cell
&
ce
=
lattice
->
getCell
(
c
);
if
(
!
ce
.
isEmpty
()
&&
!
ce
.
isFrozen
()){
uint
idx
=
ce
.
stateVectorIndex
();
m_delta
[
c
]
=
MAX
(
0
,
coupling
.
g_func
(
y
[
idx
]));
}
}
#endif
}
template
<
class
Coupling
>
void
CollierSignaling
<
Coupling
>::
fromState
(
const
state_type
&
state
){
int
n
=
lattice
->
getNbCells
();
for
(
int
c
=
0
;
c
<
n
;
c
++
){
Cell
&
ce
=
lattice
->
getCell
(
c
);
uint
idx
=
ce
.
stateVectorIndex
();
m_notch
[
c
]
=
state
[
idx
];
#ifndef ADIABATIC
m_delta
[
c
]
=
state
[
idx
+
1
];
#endif
}
}
template
<
class
Coupling
>
uint
CollierSignaling
<
Coupling
>::
cellStateDim
()
const
{
#ifdef ADIABATIC
return
1
;
#else
return
2
;
#endif
}
template
<
class
Coupling
>
float
CollierSignaling
<
Coupling
>::
getWeight
(
Boundary
*
b
,
const
state_type
&
y
){
double
diff
=
fabs
(
y
[
b
->
cell
()
->
stateVectorIndex
()]
-
y
[
b
->
oppositeCell
()
->
stateVectorIndex
()]);
//double diff = MAX(y[b->cell()->stateVectorIndex()],y[b->oppositeCell()->stateVectorIndex()]);
double
h
=
2
;
// was 2
double
aa
=
pow
(
h
,
3
);
double
bb
=
pow
(
diff
,
3
);
return
0.1
+
aa
/
(
aa
+
bb
);
// if(isHC(b->cell())){
// return isHC(b->oppositeCell())?50.0:weightfunc(m_signal[b->oppositeCell()->id()]) ;
// }
// else{
// return isSC(b->oppositeCell())?30.0:weightfunc(m_signal[b->cell()->id()]);
// }
//return 1.0;
}
template
<
class
Coupling
>
void
CollierSignaling
<
Coupling
>::
initState
(
state_type
&
state
,
double
*
signal
){
uint
nbcells
=
lattice
->
getNbCells
();
// m_notch.resize(nbcells);// already initialized
// m_delta.resize(nbcells);
for
(
uint
c
=
0
;
c
<
nbcells
;
c
++
){
Cell
&
ce
=
lattice
->
getCell
(
c
);
m_notch
[
c
]
=
state
[
ce
.
stateVectorIndex
()]
=
signal
[
2
*
c
+
1
];
m_delta
[
c
]
=
signal
[
2
*
c
];
#ifndef ADIABATIC
state
[
ce
.
stateVectorIndex
()
+
1
]
=
signal
[
2
*
c
];
#endif
// cout<<"init "<<c<< " " <<m_notch[c]<<" "<<m_delta[c]<<" "<<ce.stateVectorIndex()<<endl;
}
}
template
<
class
Coupling
>
void
CollierSignaling
<
Coupling
>::
updateState
(
state_type
&
state
){
uint
nbcells
=
lattice
->
getNbCells
();
// m_notch.resize(nbcells);// already initialized
// m_delta.resize(nbcells);
for
(
uint
c
=
0
;
c
<
nbcells
;
c
++
){
Cell
&
ce
=
lattice
->
getCell
(
c
);
if
(
ce
.
isEmpty
()){
m_notch
[
c
]
=
m_delta
[
c
]
=
state
[
ce
.
stateVectorIndex
()]
=-
1.0
;
#ifndef ADIABATIC
state
[
ce
.
stateVectorIndex
()
+
1
]
=-
1.0
;
#endif
}
}
}
// first delta and then notch
template
<
class
Coupling
>
void
CollierSignaling
<
Coupling
>::
fillVector
(
double
*
v
,
int
i
)
const
{
uint
cnt
=
0
;
uint
nbcells
=
lattice
->
getNbCells
();
for
(
uint
i
=
0
;
i
<
nbcells
;
i
++
){
v
[
cnt
++
]
=
m_delta
[
i
];
v
[
cnt
++
]
=
m_notch
[
i
];
}
}
template
<
class
Coupling
>
double
CollierSignaling
<
Coupling
>::
delaminationScore
(
Boundary
*
bo
){
double
no1
=
m_notch
[
bo
->
cell
()
->
id
()];
double
no2
=
m_notch
[
bo
->
oppositeCell
()
->
id
()];
double
diff
=
fabs
(
no1
-
no2
);
return
1.0
/
diff
;
}
#endif
#ifdef TEST_COLLIER
int
main
(
int
argc
,
char
*
argv
[]){
HillCauchyCoupling
sig
;
sig
.
init
();
double
x
=
0.6
;
cout
<<
x
<<
":"
<<
sig
.
f_func
(
x
)
<<
endl
;
return
0
;
}
#endif
Event Timeline
Log In to Comment