Choral
cardioModel_mod.f90
Go to the documentation of this file.
1 
10 
12 
15  use real_type, only: rp
16  use basic_tools
17  use funclib
18  use abstract_interfaces, only: r3tor3
19  use algebra_lin
20  use fespace_mod
21  use quadmesh_mod
22  use diffusion
23  use csr_mod
24 
25  implicit none
26  private
27 
28 
29  ! %----------------------------------------%
30  ! | |
31  ! | PUBLIC DATA |
32  ! | |
33  ! %----------------------------------------%
34  public :: cardiomodel
35 
36  public :: print
37  public :: monodomain_assemble
38  public :: scale_conductivities
39 
40 
41  ! %----------------------------------------%
42  ! | |
43  ! | DERIVED TYPE |
44  ! | |
45  ! %----------------------------------------%
47  type :: cardiomodel
48 
49  !! data provided by the constructor
50  !!
52  procedure(r3tor3) , pointer, nopass :: fibre => vector_field_e_x
53 
55  real(RP) :: am
56 
57  !! Data defined by the code
58  !!
60  real(RP) :: gil
62  real(RP) :: git
64  real(RP) :: gel
66  real(RP) :: get
67 
68  end type cardiomodel
69 
70 
71  ! %----------------------------------------%
72  ! | |
73  ! | GENERIC INTERFACES |
74  ! | |
75  ! %----------------------------------------%
77  interface print
78  module procedure cardiomodel_print
79  end interface print
80 
82  interface cardiomodel
83  module procedure cardiomodel_create
84  end interface cardiomodel
85 
86 contains
87 
98  function cardiomodel_create(fibre, Am, conduc_type) result(cm)
99 
100  type(cardiomodel) :: cm
101  procedure(r3tor3) :: fibre
102  real(RP) , intent(in) :: am
103  integer , intent(in) :: conduc_type
104 
105  if (choral_verb>0) write(*,*) &
106  & "cardioModel_mod : create"
107 
108  cm%fibre => fibre
109 
110  cm%Am = am
111 
112  call set_conductivities(cm%GIL, &
113  &cm%GIT, cm%GEL, cm%GET, conduc_type)
114 
115  end function cardiomodel_create
116 
117 
120  subroutine set_conductivities(gil, git, gel, get, type)
121  real(RP), intent(out) :: gil, git,gel,get
122  integer , intent(in) :: type
123 
124  if (choral_verb>1) write(*,"(A36,$)") &
125  & " cardioModel_mod : set_conduct. = "
126  select case(type)
127  case(le_guyader)
128  if (choral_verb>1) write(*,*) "LE_GUYADER"
129  gil = 1.741_rp
130  git = 0.1934_rp
131  gel = 3.906_rp
132  get = 1.970_rp
133 
134  case(le_guyader_2)
135  if (choral_verb>1) write(*,*) "LE_GUYADER_2"
136  gil = 1.741_rp
137  git = 0.475_rp
138  gel = 3.906_rp
139  get = 1.970_rp
140 
141  case(cond_iso_1)
142  if (choral_verb>1) write(*,*) "COND_ISO_1"
143  gil = 1._rp
144  git = 1._rp
145  gel = 1._rp
146  get = 1._rp
147 
148  case default
149  call quit("cardioModel_mod: set_conductivities: 1")
150 
151  end select
152 
153  end subroutine set_conductivities
154 
155 
158  subroutine scale_conductivities(cm, sc)
159  type(cardiomodel), intent(inout) :: cm
160  real(RP) , intent(in) :: sc
161 
162  cm%GIL = cm%GIL * sc
163  cm%GIT = cm%GIT * sc
164  cm%GEL = cm%GEL * sc
165  cm%GET = cm%GET * sc
166 
167  end subroutine scale_conductivities
168 
169 
172  subroutine cardiomodel_print(cm)
173  type(cardiomodel), intent(in) :: cm
174 
175  write(*,*)"cardioModel_mod : print "
176  write(*,*)" Am =", cm%am
177  write(*,*)" Intra conduc. GIL =", cm%GIL
178  write(*,*)" GIT =", cm%GIT
179  write(*,*)" Extra conduc. GEL =", cm%GEL
180  write(*,*)" GET =", cm%GET
181 
182  end subroutine cardiomodel_print
183 
184 
197  subroutine monodomain_assemble(M, S, cm, X_h, qd_M, qd_S)
198  type(csr) , intent(out):: M, S
199  class(cardiomodel), intent(in) :: cm
200  type(fespace) , intent(in) :: X_h
201  integer , intent(in) :: qd_M, qd_S
202 
203  type(quadmesh) :: qdm
204  real(RP) :: gl, gt
205 
206  if (choral_verb>0) write(*,*) &
207  & "cardioModel_mod : monodomain"
208 
209  ! monodomain conductivities
210  gl = 2._rp / (1._rp/cm%GIL + 1._rp/cm%GEL )
211  gt = 2._rp / (1._rp/cm%GIT + 1._rp/cm%GET )
212  if (choral_verb>1) write(*,*) &
213  & " conduc. long/trans =", &
214  & real(gl,SP), real(gt,SP)
215 
216  ! mass matrix
217  qdm = quadmesh(x_h%mesh)
218  call set(qdm, qd_m)
219  call assemble(qdm)
220  call diffusion_massmat(m, one_r3, x_h, qdm)
221 
222  ! stiffness matrix
223  qdm = quadmesh(x_h%mesh)
224  call set(qdm, qd_s)
225  call assemble(qdm)
226  call diffusion_stiffmat(s, b, x_h, qdm)
227  call scale(s, 1._rp/cm%Am)
228 
229  contains
230 
231  function b(x, w1, w2)
232  real(RP) :: b
233  real(RP), dimension(3), intent(in) :: x, w1, w2
234 
235  real(RP), dimension(3) :: fd
236 
237  fd = cm%fibre(x)
238  b = usigmav(w1, w2, fd, gl, gt)
239 
240  end function b
241 
242  end subroutine monodomain_assemble
243 
244 
245 
255  function usigmav(u, v, fd, gl, gt) result(r)
256  real(RP) :: r
257  real(RP), dimension(3), intent(in) :: u, v, fd
258  real(RP) , intent(in) :: gl, gt
259 
260  real(RP), dimension(3) :: e1, e2, e3
261  real(RP) :: nrm
262 
263  nrm = sqrt(sum(fd*fd))
264  e1 = fd / nrm
265 
266  e2(1) = e1(3) - e1(2)
267  e2(2) = e1(1) - e1(3)
268  e2(3) = e1(2) - e1(1)
269 
270  nrm = sqrt(sum(e2*e2))
271  e2 = e2/nrm
272 
273  e3 = crossprod_3d(e1,e2)
274 
275  r = sum(u*e1) * sum(v*e1) * gl
276  r = r + sum(u*e2) * sum(v*e2) * gt
277  r = r + sum(u*e3) * sum(v*e3) * gt
278 
279  end function usigmav
280 
281 
282 end module cardiomodel_mod
DEFINITION OF ABSTRACT_INTERFACES FOR THE LIBRARY CHORAL
scale the matrix entries
Definition: csr_mod.F90:150
LINEAR ALGEBRA TOOLS
Definition: algebra_lin.F90:9
BASIC TOOLS
Definition: basic_tools.F90:8
integer, parameter, public rp
real(kind=RP) = real precision in the code REAL_TOL = epsilon to test real equality ...
Definition: real_type.F90:90
Derived type for tissue properties in cardiac electrophysiology.
integer, parameter le_guyader_2
subroutine, public quit(message)
Stop program execution, display an error messahe.
real(rp) function, public one_r3(x)
The function .
Definition: funcLib.f90:67
DERIVED TYPE cardioModel
PRE-DEFINED NUMERIC FUNCTIONS
Definition: funcLib.f90:25
subroutine, public diffusion_massmat(mass, a, X_h, qdm, dofToDof)
Assemble the mass matrix of the bilinear product:
Definition: diffusion.F90:168
real(rp) function, dimension(3), public vector_field_e_x(x)
constant vector field
Definition: funcLib.f90:56
REAL NUMBERS PRECISION IN CHORAL: selects simple/double/quad
Definition: real_type.F90:33
DERIVED TYPE quadMesh: integration methods on meshes
DERIVED TYPE feSpace: finite element spaces
Definition: feSpace_mod.F90:58
MODULE FOR DIFFUSION PROBLEMS
Definition: diffusion.F90:24
type(cardiomodel) function cardiomodel_create(fibre, Am, conduc_type)
constructor for cardioModel type
CHORAL CONSTANTS
real(rp) function usigmav(u, v, fd, gl, gt)
scalar product u^T Sigma v
Derived type for sparse matrices in CSR format.
Definition: csr_mod.F90:62
print a short description
integer, parameter cond_iso_1
Derived type feSpace: finite element space on a mesh.
integer choral_verb
Verbosity level.
subroutine, public scale_conductivities(cm, sc)
rescale the conductivities
subroutine, public diffusion_stiffmat(stiff, b, X_h, qdm, dofToDof)
Assemble the stiffness matrix of the bilinear product:
Definition: diffusion.F90:366
integer, parameter le_guyader
DEFINITION OF GLOBAL VARIABLES FOR THE LIBRARY CHORAL
The type quadMesh defines integration methods on meshes.
real(rp) function, dimension(3), public crossprod_3d(u, v)
Cross product in R**3.
subroutine cardiomodel_print(cm)
print cardio model
DERIVED TYPE csr for sparse matrices
Definition: csr_mod.F90:13
subroutine set_conductivities(gil, git, gel, get, type)
set conductivities
subroutine, public monodomain_assemble(M, S, cm, X_h, qd_M, qd_S)
Assemble the monodomain diffusion matrices.
real(rp) function b(x, w1, w2)