Academia.eduAcademia.edu

Figure 5 – uploaded by Martin Bishop

See full PDFdownloadDownload figure
Fig. 4. Data from 20 ms simulations on a 2cm 1D domain using our method with and without the Noble-form LR91 cell model modification. Note th limited accuracy caused by the discontinuity in standard LR91. At= 0.001 ms. L*(L) is the maximum-in-time of the L?-norm of the error in space. Th errors are against a quartic reference solution with h = 10-4 cm.

Figure 4 Data from 20 ms simulations on a 2cm 1D domain using our method with and without the Noble-form LR91 cell model modification. Note th limited accuracy caused by the discontinuity in standard LR91. At= 0.001 ms. L*(L) is the maximum-in-time of the L?-norm of the error in space. Th errors are against a quartic reference solution with h = 10-4 cm.

Related Figures (17)

Fig. 1. Three of the ten hierarchical basis function pieces required in 2D on the reference element for a cubic finite element approximation of the solution. Note that when mapped to the real mesh from the reference element, each of these is just part of one of the basis functions that is actually supported on multiple elements. We begin with a necessary lemma, supposing that the true solution to system (2) has at least k derivatives.
Fig. 1. Three of the ten hierarchical basis function pieces required in 2D on the reference element for a cubic finite element approximation of the solution. Note that when mapped to the real mesh from the reference element, each of these is just part of one of the basis functions that is actually supported on multiple elements. We begin with a necessary lemma, supposing that the true solution to system (2) has at least k derivatives.
Fig. 2. 1D simulation results demonstrating that we never do better than quadratic convergence in the L~(L”) norm when the approximation to w is linear- only, even when the finite element space S),, is of high enough order to allow for better convergence. The points are the measured errors, the dashed line shows a theoretical quadratic convergence gradient. The simulation was 1D and the errors are measured against a reference with mesh spacing h=0.0001 cm, whereas the test simulations use meshes with six different spacings ranging from h = 0.01 cm to h = 0.001 cm.
Fig. 2. 1D simulation results demonstrating that we never do better than quadratic convergence in the L~(L”) norm when the approximation to w is linear- only, even when the finite element space S),, is of high enough order to allow for better convergence. The points are the measured errors, the dashed line shows a theoretical quadratic convergence gradient. The simulation was 1D and the errors are measured against a reference with mesh spacing h=0.0001 cm, whereas the test simulations use meshes with six different spacings ranging from h = 0.01 cm to h = 0.001 cm.
the L7(Q) projection into S,, in each of its m+1 components, with m the number of components of w, and also the (m+1) x (m+1) matrix-like operator The importance of p can be seen via an a priori error estimate (see [41]). We examine this in what follows, where we work with the vectorised formulation of system (2), treating the spatial discretisation of the non-diffusing state variables in the finite element framework. To this end we introduce the inner product ([a12]", [b1b2]"), := (a1, b1) + 32, (a2, bz;) and asso- ciated norm || - ||-, where (a,b) is the standard L? inner product, a;, by € L?(Q) and dy, by € (L7(Q))™. Define
the L7(Q) projection into S,, in each of its m+1 components, with m the number of components of w, and also the (m+1) x (m+1) matrix-like operator The importance of p can be seen via an a priori error estimate (see [41]). We examine this in what follows, where we work with the vectorised formulation of system (2), treating the spatial discretisation of the non-diffusing state variables in the finite element framework. To this end we introduce the inner product ([a12]", [b1b2]"), := (a1, b1) + 32, (a2, bz;) and asso- ciated norm || - ||-, where (a,b) is the standard L? inner product, a;, by € L?(Q) and dy, by € (L7(Q))™. Define
Fig. 3. Plot comparing the action potential in an isolated cell model using the standard Luo-Rudy 1991 formulation and the Noble-form modification. Left Noble-form modified LR91 action potential in an isolated cell. On this scale, differences caused by the modification would be hardly noticeable. Right original LR91 transmembrane potential subtracted from Noble-form LR91 transmembrane potential. Note that the scales on the y-axes differ. We stimulated a 2 cm 1D domain at one end with a ramp stimulus using a time-step of At = 0.001 ms and simulated the first 20 ms of activation. The errors in two different norms are presented in Fig. 4(b) and Fig. 5(a) and use as a reference a quartic solution generated with h = 0.0001 cm. Fig. 5(a) shows the error measured using the L?-in-time norm of the #!(Q) norm in space, for which we expect O(h?) convergence gradients [15]. Note the agreement of Figs. 4(b) and 5(a) with the
Fig. 3. Plot comparing the action potential in an isolated cell model using the standard Luo-Rudy 1991 formulation and the Noble-form modification. Left Noble-form modified LR91 action potential in an isolated cell. On this scale, differences caused by the modification would be hardly noticeable. Right original LR91 transmembrane potential subtracted from Noble-form LR91 transmembrane potential. Note that the scales on the y-axes differ. We stimulated a 2 cm 1D domain at one end with a ramp stimulus using a time-step of At = 0.001 ms and simulated the first 20 ms of activation. The errors in two different norms are presented in Fig. 4(b) and Fig. 5(a) and use as a reference a quartic solution generated with h = 0.0001 cm. Fig. 5(a) shows the error measured using the L?-in-time norm of the #!(Q) norm in space, for which we expect O(h?) convergence gradients [15]. Note the agreement of Figs. 4(b) and 5(a) with the
Fig. 5. Noble-form LR91 with our method. Fig. 5(a) shows the L?(H') norm of the error; this is the L?-in-time norm of the Sobolev || - ||1;,2 norm of the error in space. Fig. 5(b) shows the same data as Fig. 4(b), but with the exponential convergence in p highlighted instead. Data from 20 ms simulations on a 2 cm 1D domain. At = 0.001 ms. The errors are against a quartic reference solution with h = 10-4 cm. The values of h given are in cm.
Fig. 5. Noble-form LR91 with our method. Fig. 5(a) shows the L?(H') norm of the error; this is the L?-in-time norm of the Sobolev || - ||1;,2 norm of the error in space. Fig. 5(b) shows the same data as Fig. 4(b), but with the exponential convergence in p highlighted instead. Data from 20 ms simulations on a 2 cm 1D domain. At = 0.001 ms. The errors are against a quartic reference solution with h = 10-4 cm. The values of h given are in cm.
Six cm 1D domain simulation with a ramp stimulus at one end and At = 0.01 ms. The activation time is the time until the transmembrane potential at the nod at the opposite end of the domain (at 6 cm) passes up through 0 mV. CVs are computed using the nodes at 2 and 4cm.
Six cm 1D domain simulation with a ramp stimulus at one end and At = 0.01 ms. The activation time is the time until the transmembrane potential at the nod at the opposite end of the domain (at 6 cm) passes up through 0 mV. CVs are computed using the nodes at 2 and 4cm.
Fig. 6. Wavefront location at t= 12 ms with p= 1-4, demonstrating the scheme with isotropic tissue and with fibres aligned with the x-axis having anisotropic conductivity (1S m7! along the fibres and 0.2 S m~! perpendicular to them). In both cases, the larger p is, the less distance the wavefront has propagated. Note that p = 2-4 are indistinguishable in (a), and that p = 3 and p = 4 are indistinguishable in (b). This simulation used At = 0.01 ms, the mean element diameter was 0.0113 cm and the stimulus was applied to the lower-left corner of the domain.
Fig. 6. Wavefront location at t= 12 ms with p= 1-4, demonstrating the scheme with isotropic tissue and with fibres aligned with the x-axis having anisotropic conductivity (1S m7! along the fibres and 0.2 S m~! perpendicular to them). In both cases, the larger p is, the less distance the wavefront has propagated. Note that p = 2-4 are indistinguishable in (a), and that p = 3 and p = 4 are indistinguishable in (b). This simulation used At = 0.01 ms, the mean element diameter was 0.0113 cm and the stimulus was applied to the lower-left corner of the domain.
Activation times, conduction velocities, percentage conduction velocity errors and linear system solve times for a variety of 2D meshes of Q = [0,1] x [0,1] cm. The first timings presented are for a single time-step only and use preconditioned conjugate gradients (PCG) with an incomplete LU (ILU) decomposition of the portion of the system matrix corresponding to the linear basis functions as a preconditioner, and the second times all code on a time-step iteration (ILU PCG + cell updates, etc.). The stimulus was along the bottom edge of the domain and the activation time is for the node in the top-right corner. Asterisked basis degrees indicate that the simulation was run with A t = 0.001 ms instead of the usual At = 0.01 ms. The percentage errors use as a reference the p = 4 simulation on the finest mesh with the same value of At. results. The CV is an average taken over many randomly selected point pairs in the domain. Percentage errors in conduction velocity are also presented in the table, using the finest quartic simulation with the same value of At as a reference. Note the high accuracy of the quartic simulations, and the minimal variation in the percentage errors caused by varying At. Details on the degrees of freedom are in Table 5.
Activation times, conduction velocities, percentage conduction velocity errors and linear system solve times for a variety of 2D meshes of Q = [0,1] x [0,1] cm. The first timings presented are for a single time-step only and use preconditioned conjugate gradients (PCG) with an incomplete LU (ILU) decomposition of the portion of the system matrix corresponding to the linear basis functions as a preconditioner, and the second times all code on a time-step iteration (ILU PCG + cell updates, etc.). The stimulus was along the bottom edge of the domain and the activation time is for the node in the top-right corner. Asterisked basis degrees indicate that the simulation was run with A t = 0.001 ms instead of the usual At = 0.01 ms. The percentage errors use as a reference the p = 4 simulation on the finest mesh with the same value of At. results. The CV is an average taken over many randomly selected point pairs in the domain. Percentage errors in conduction velocity are also presented in the table, using the finest quartic simulation with the same value of At as a reference. Note the high accuracy of the quartic simulations, and the minimal variation in the percentage errors caused by varying At. Details on the degrees of freedom are in Table 5.
Activation times in the 1 cm by 1 cm domain with fibre orientation and four measurement points shown in Fig. 7(a) using various h and p values fine structure that must be properly captured using small elements. This vector field was generated using the Laplace- Dirichlet approach and the time-step used was At = 0.01 ms. The activation times at the three points shown in Fig. 7(b are given in Table 4. Our results show that the second coarsest mesh with p = 3 can half the worst percentage error in acti- vation time using a linear system which can be solved twice as quickly when compared to p = 1 on the finest mesh, or with p =4, we can reduce the error to one-tenth that of p = 1 on the finest mesh with a linear system which takes slightly longer tc solve. See Fig. 8 for some wavefront locations using various meshes and values of p; note how poor the wavefront location can be when p = 1 even on highly refined meshes, and how on the second coarsest mesh with p = 3 the results are better thar the computationally twice-as-demanding p = 1 on the finest mesh. Table 3
Activation times in the 1 cm by 1 cm domain with fibre orientation and four measurement points shown in Fig. 7(a) using various h and p values fine structure that must be properly captured using small elements. This vector field was generated using the Laplace- Dirichlet approach and the time-step used was At = 0.01 ms. The activation times at the three points shown in Fig. 7(b are given in Table 4. Our results show that the second coarsest mesh with p = 3 can half the worst percentage error in acti- vation time using a linear system which can be solved twice as quickly when compared to p = 1 on the finest mesh, or with p =4, we can reduce the error to one-tenth that of p = 1 on the finest mesh with a linear system which takes slightly longer tc solve. See Fig. 8 for some wavefront locations using various meshes and values of p; note how poor the wavefront location can be when p = 1 even on highly refined meshes, and how on the second coarsest mesh with p = 3 the results are better thar the computationally twice-as-demanding p = 1 on the finest mesh. Table 3
Fig. 7. Fibre orientations used for the inhomogeneous simulations. The marked points are those at which activation time is measured. For Fig. 7(a), results are displayed in Table 3 and stimulation was a ramp in the bottom-left corner. For Fig. 7(b), results are displayed in in Table 4 and stimulation was a ramp along the bottom of the domain.
Fig. 7. Fibre orientations used for the inhomogeneous simulations. The marked points are those at which activation time is measured. For Fig. 7(a), results are displayed in Table 3 and stimulation was a ramp in the bottom-left corner. For Fig. 7(b), results are displayed in in Table 4 and stimulation was a ramp along the bottom of the domain.
Activation times in the 1 cm x 1 cm domain with holes, fibre orientation and three measurement points shown in Fig. 7(b) using various h and p values. The meshes are identified by maximum element diameter because the mean element size would convey little information due to the very small elements near the holes. The timings presented are for the linear system only (ILU PCG) and for the whole-timestep (ILU PCG + cell updates, etc.). More information about the meshes used is presented in Table 5, and some plots of the wavefront at 17 ms are shown in Fig. 8.
Activation times in the 1 cm x 1 cm domain with holes, fibre orientation and three measurement points shown in Fig. 7(b) using various h and p values. The meshes are identified by maximum element diameter because the mean element size would convey little information due to the very small elements near the holes. The timings presented are for the linear system only (ILU PCG) and for the whole-timestep (ILU PCG + cell updates, etc.). More information about the meshes used is presented in Table 5, and some plots of the wavefront at 17 ms are shown in Fig. 8.
Table showing the degrees of freedom used throughout this study.
Table showing the degrees of freedom used throughout this study.
Fig. 8. Wavefront location at 17 ms using various meshes and degrees. (a) Transmembrane potential distribution on the second coarsest mesh with p = 4. Note that each triangle represents a single element of the mesh, and that this is independent of p. (b) Wavefront (u = 0) contours on the second coarsest mesh with various degrees compared to p = 1 on the finest mesh. (c) Wavefront contours on the second coarsest mesh with various degrees compared to p=1 on the finest mesh. (d) Detail from Fig. 8(c). See the corresponding activation time data in Table 4.
Fig. 8. Wavefront location at 17 ms using various meshes and degrees. (a) Transmembrane potential distribution on the second coarsest mesh with p = 4. Note that each triangle represents a single element of the mesh, and that this is independent of p. (b) Wavefront (u = 0) contours on the second coarsest mesh with various degrees compared to p = 1 on the finest mesh. (c) Wavefront contours on the second coarsest mesh with various degrees compared to p=1 on the finest mesh. (d) Detail from Fig. 8(c). See the corresponding activation time data in Table 4.
(b) Two holes mesh, second coarsest, p=1-4 and finest mesh p=1.
(b) Two holes mesh, second coarsest, p=1-4 and finest mesh p=1.
Simulations on the mesh with holes in the context of the prescribed error tolerances for percentage activation time error at which they would be acceptable. [he average single-iteration linear system solve time ratio (time for simulation)/(time for finest linear simulation) is given. Finest mesh p = 4 solution taken to be the “true” solution. Errors are worst-case over the three test nodes used. Full data in Table 4. imulations using the very simple FitzHugh-Nagumo model [48] instead of LR91 with p = 4 do not display any problems, either do p = 6 simulations with LR91 for p = 1-4. Regardless of this, the efficiency remains compelling, and we note that the time-integrated Sobolev norm spatial error Fig. 5(a)) does not appear to be so strongly affected.
Simulations on the mesh with holes in the context of the prescribed error tolerances for percentage activation time error at which they would be acceptable. [he average single-iteration linear system solve time ratio (time for simulation)/(time for finest linear simulation) is given. Finest mesh p = 4 solution taken to be the “true” solution. Errors are worst-case over the three test nodes used. Full data in Table 4. imulations using the very simple FitzHugh-Nagumo model [48] instead of LR91 with p = 4 do not display any problems, either do p = 6 simulations with LR91 for p = 1-4. Regardless of this, the efficiency remains compelling, and we note that the time-integrated Sobolev norm spatial error Fig. 5(a)) does not appear to be so strongly affected.
Related topics:
EngineeringComputational PhysicsFinite Element MethodsMonodomain Simulation
580 California St., Suite 400
San Francisco, CA, 94104
© 2025 Academia. All rights reserved