OpenCMISS Demonstrations

The purpose of this effort is to document some of the OpenCMISS examples as demonstrations of some of the features available in the OpenCMISS platform. The actual examples often combine multiple aspects of several important features, and as such individual examples may be found under several of the top-level feature categories introduced below.

Runtime definition of mathematical models

OpenCMISS makes use of the CellML format to enable users to define computational fields in their model(s) to be defined based on mathematical models that are specified at the time at which a simulation is executed. Such fields may be used in many different applications and the demonstrations below highlight some of the current examples.

OpenCMISS CellML Examples

These demonstrations illustrate the basic creation of computational fields which are defined using CellML models in each of the primary language bindings provided by the OpenCMISS library.

Fortran introduction to CellML in OpenCMISS

This example provides an entry level demonstration of creating fields in OpenCMISS to be defined using CellML models. This demonstration uses Fortran, see Python introduction to CellML in OpenCMISS for the corresponding Python example. The source code for this example is available here: src/FortranExample.f90.

Following the usual OpenCMISS practices, you first need to declare the objects to be used in you application. For CellML, we would normally declare:

1
2
  TYPE(CMISSCellMLType) :: CellML
  TYPE(CMISSFieldType) :: CellMLModelsField,CellMLStateField,CellMLIntermediateField,CellMLParametersField

which declares a single CellML environment that we will use and the fields that we will use with CellML. We also fetch the URL of the CellML model we will load from the command line arguments, if it is provided. Otherwise we fall back on the default Noble 1998 model.

1
2
3
4
5
6
7
8
  ! we want to get the CellML model from the command line
  NUMBER_OF_ARGUMENTS = COMMAND_ARGUMENT_COUNT()
  IF(NUMBER_OF_ARGUMENTS >= 1) THEN
    CALL GET_COMMAND_ARGUMENT(1,ModelUrl,ARGUMENT_LENGTH,STATUS)
    IF(STATUS>0) WRITE(*,'(">>ERROR: Error for command argument 1.")')
  ELSE
    ModelUrl="n98.xml"
  ENDIF

We then create our CellML environment and import the model defined above. Having imported the model we are able to flag the parameters of interest, namely the stimulus current and the sodium channel conductance, as known variables. We also flag all the membrane currents as wanted variables in order to have access to their values from OpenCMISS.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
  !Create the CellML environment
  CALL CMISSCellML_Initialise(CellML,Err)
  CALL CMISSCellML_CreateStart(CellMLUserNumber,Region,CellML,Err)
  !Import the specified model into the CellML environment
  CALL CMISSCellML_ModelImport(CellML,ModelUrl,modelIndex,Err)
  ! Now we have imported all the models we are able to specify which variables from the model we want:
  !   - to set from this side
  CALL CMISSCellML_VariableSetAsKnown(CellML,modelIndex,"fast_sodium_current/g_Na ",Err)
  CALL CMISSCellML_VariableSetAsKnown(CellML,modelIndex,"membrane/IStim",Err)
  !   - to get from the CellML side (state variables are wanted by default and can not be changed)
  CALL CMISSCellML_VariableSetAsWanted(CellML,modelIndex,"membrane/i_K1",Err)
  CALL CMISSCellML_VariableSetAsWanted(CellML,modelIndex,"membrane/i_to",Err)
  CALL CMISSCellML_VariableSetAsWanted(CellML,modelIndex,"membrane/i_K",Err)
  CALL CMISSCellML_VariableSetAsWanted(CellML,modelIndex,"membrane/i_K_ATP",Err)
  CALL CMISSCellML_VariableSetAsWanted(CellML,modelIndex,"membrane/i_Ca_L_K",Err)
  CALL CMISSCellML_VariableSetAsWanted(CellML,modelIndex,"membrane/i_b_K",Err)
  CALL CMISSCellML_VariableSetAsWanted(CellML,modelIndex,"membrane/i_NaK",Err)
  CALL CMISSCellML_VariableSetAsWanted(CellML,modelIndex,"membrane/i_Na",Err)
  CALL CMISSCellML_VariableSetAsWanted(CellML,modelIndex,"membrane/i_b_Na",Err)
  CALL CMISSCellML_VariableSetAsWanted(CellML,modelIndex,"membrane/i_Ca_L_Na",Err)
  CALL CMISSCellML_VariableSetAsWanted(CellML,modelIndex,"membrane/i_NaCa",Err)
  !   - and override constant parameters without needing to set up fields
  !> \todo Need to allow parameter values to be overridden for the case when user has non-spatially varying parameter value.
!  CALL CMISSDiagnosticsSetOff(Err)
  !Finish the CellML environment
  CALL CMISSCellML_CreateFinish(CellML,Err)

Having flagged the desired variables, we finish the creation of the CellML environment. This triggers the instantion of the CellML model into a computable black box for use by solvers in OpenCMISS control loops.

We are then able to specify the actual field mapping, in this case assuming that we are going to perform an electrophysiology simulation with the CellML model so needing to map the membrane potential variable from the model to the OpenCMISS dependent field.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
  !Start the creation of CellML <--> OpenCMISS field maps
  CALL CMISSCellML_FieldMapsCreateStart(CellML,Err)
  !Now we can set up the field variable component <--> CellML model variable mappings.
  !Map Vm - field to CellML
  CALL CMISSCellML_CreateFieldToCellMLMap(CellML,DependentField,CMISS_FIELD_U_VARIABLE_TYPE,1,CMISS_FIELD_VALUES_SET_TYPE, &
    & modelIndex,"membrane/V",CMISS_FIELD_VALUES_SET_TYPE,Err)
  ! - and CellML to field
  CALL CMISSCellML_CreateCellMLToFieldMap(CellML,modelIndex,"membrane/V",CMISS_FIELD_VALUES_SET_TYPE, &
    & DependentField,CMISS_FIELD_U_VARIABLE_TYPE,1,CMISS_FIELD_VALUES_SET_TYPE,Err)
  !Finish the creation of CellML <--> OpenCMISS field maps
  CALL CMISSCellML_FieldMapsCreateFinish(CellML,Err)

Following the maps we create the CellML fields:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
  !Start the creation of the CellML models field
  CALL CMISSField_Initialise(CellMLModelsField,Err)
  CALL CMISSCellML_ModelsFieldCreateStart(CellML,CellMLModelsFieldUserNumber,CellMLModelsField,Err)
  !Finish the creation of the CellML models field
  CALL CMISSCellML_ModelsFieldCreateFinish(CellML,Err)

  !Start the creation of the CellML state field
  CALL CMISSField_Initialise(CellMLStateField,Err)
  CALL CMISSCellML_StateFieldCreateStart(CellML,CellMLStateFieldUserNumber,CellMLStateField,Err)
  !Finish the creation of the CellML state field
  CALL CMISSCellML_StateFieldCreateFinish(CellML,Err)

  !Start the creation of the CellML intermediate field
  CALL CMISSField_Initialise(CellMLIntermediateField,Err)
  CALL CMISSCellML_IntermediateFieldCreateStart(CellML,CellMLIntermediateFieldUserNumber,CellMLIntermediateField,Err)
  !Finish the creation of the CellML intermediate field
  CALL CMISSCellML_IntermediateFieldCreateFinish(CellML,Err)
  
  !Start the creation of CellML parameters field
  CALL CMISSField_Initialise(CellMLParametersField,Err)
  CALL CMISSCellML_ParametersFieldCreateStart(CellML,CellMLParametersFieldUserNumber,CellMLParametersField,Err)
  !Finish the creation of CellML parameters
  CALL CMISSCellML_ParametersFieldCreateFinish(CellML,Err)

We then define the spatial variation for the stimulus current variable flagged as known above. When setting up this spatial distribution we have to take into account the fact that this example application may be run in parallel and therefore only set values relevant to each specific computational node.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
  CALL CMISSCellML_FieldComponentGet(CellML,modelIndex,CMISS_CELLML_PARAMETERS_FIELD,"membrane/IStim",stimcomponent,Err)
  !Set the Stimulus at half the bottom nodes
  DO node_idx=1,NUMBER_OF_ELEMENTS/2
    CALL CMISSDecomposition_NodeDomainGet(Decomposition,node_idx,1,nodeDomain,Err)
    IF(nodeDomain==ComputationalNodeNumber) THEN
      CALL CMISSField_ParameterSetUpdateNode(CellMLParametersField,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,1, &
        & node_idx, &
        & stimcomponent,STIM_VALUE,Err)
    ENDIF
  ENDDO

And finally, in a similar manner as the electrical stimulus above, we define the spatial distribution of the sodium channel conductance variable also flagged as known above.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
  !Set up the g_Na gradient
  CALL CMISSCellML_FieldComponentGet(CellML,modelIndex,CMISS_CELLML_PARAMETERS_FIELD,"fast_sodium_current/g_Na", &
    & gNacomponent,Err)
  !Loop over the nodes
  DO node_idx=1,(NUMBER_OF_ELEMENTS+1)*(NUMBER_OF_ELEMENTS+1)
    CALL CMISSField_ParameterSetGetNode(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,1,node_idx,1, &
      & X,Err)
    CALL CMISSField_ParameterSetGetNode(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,1,node_idx,2, &
      & Y,Err)
    DISTANCE=SQRT(X**2+Y**2)/SQRT(2.0_CMISSDP)
    gNa_VALUE=2.0_CMISSDP*(DISTANCE+0.5_CMISSDP)*385.5e-3_CMISSDP
    CALL CMISSField_ParameterSetUpdateNode(CellMLParametersField,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,1, &
      & node_idx, &
      & gNacomponent,gNa_VALUE,Err)
  ENDDO

Python introduction to CellML in OpenCMISS

This example provides an entry level demonstration of creating fields in OpenCMISS to be defined using CellML models. This demonstration uses Python, see Fortran introduction to CellML in OpenCMISS for the correpsponding Fortran example.

Axial extension in a homogeneous pipe

Section author: Jagir Hussan <https://unidirectory.auckland.ac.nz/profile/rjag008>, David Nickerson <http://about.me/david.nickerson/>

This example demonstrates the application of an mechanical extension to a homogeneous tissue cylinder model along the axis of the cylinder. The Mooney-Rivlin constitutive law is used via CellML.

This is a Python example. In order to exectute this example, you will need to have the OpenCMISS Python bindings available. This Cmgui command file can be used to visualise the simulation output. The geometric finite element mesh used in this example was created using PyZinc with the script quadcylinderthreematerials.py.

The example is encoded in HomogeneousPipeAxialExtension.py. Below we describe some of the key aspects of this example.

Example description

In this example the The Python script exfile.py is used to load the geometry from a file in the Zinc EX format.

1
2
3
import sys, os, exfile
# Make sure $OPENCMISS_ROOT/cm/bindings/python is first in our PYTHONPATH.
sys.path.insert(1, os.path.join((os.environ['OPENCMISS_ROOT'],'cm','bindings','python')))

The exfile.py script in needs to be in the current folder or available in your Python environment for the import on line 1 to succeed. Whereas on line 3 we are specifically adding the location we expect to find the OpenCMISS Python bindings using the OPENCMISS_ROOT environment variable. Using the exfile module we are able to load the finite element model geometry from the EXREGION file: hetrogenouscylinder.exregion, as shown below:

1
2
#Load the mesh information in the form of exregion format
exregion = exfile.Exregion("hetrogenouscylinder.exregion")

Throughout the remainder of the Python script you can see the data from the now defined exregion object used to create the geometric mesh via the standard OpenCMISS methods. For example, in this code:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
# Read the geometric field 
for node_num in exregion.nodeids:
    version = 1
    derivative = 1
    coord = []
    for component in range(1, numberOfXi + 1):
        component_name = ["1", "2", "3"][component - 1]
        value = exregion.node_value("coordinates", component_name, node_num, derivative)
        geometricField.ParameterSetUpdateNodeDP(
                CMISS.FieldVariableTypes.U,
                CMISS.FieldParameterSetTypes.VALUES,
                version, derivative, node_num, component, value)
        coord.append(value)
#The cylinder has an axial length of 5 units and coord[2] contains the axial coordinate value        
    if coord[2] < 0.001:
        left_boundary_nodes.append(node_num)
    elif coord[2] > 4.999:
        right_boundary_nodes.append(node_num)
           

geometricField.ParameterSetUpdateFinish(CMISS.FieldVariableTypes.U,CMISS.FieldParameterSetTypes.VALUES)

we set the nodal coordinates and derivatives in the models geometric field. There is also code to keep track of nodes located at either end of the vessel, which will be used later when defining boundary conditions.

CellML

As mentioned above, this example uses the Mooney-Rivlin mechanical constitutive law defined in CellML. In this section we briefly describe how this is achieved, further details and general description can be found in the OpenCMISS CellML Examples. First the CellML environment must be created and our Mooney-Rivlin model imported into the environment:

1
2
3
4
5
# Create the CellML environment
CellML = CMISS.CellML()
CellML.CreateStart(CellMLUserNumber, region)
# Import a Mooney-Rivlin material law from a file
MooneyRivlinModel = CellML.ModelImport("mooney_rivlin.xml")

With the CellML model imported, we are now able to flag the variables from the model that are known (by the finite elasticity model) and those variables that are wanted (by the finite elasticity model) from the CellML model. These are the variables from the CellML model that we will want to associate with fields in the OpenCMISS (continuum) model.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
# Now we have imported the model we are able to specify which variables from the model we want to set from openCMISS
CellML.VariableSetAsKnown(MooneyRivlinModel, "equations/E11")
CellML.VariableSetAsKnown(MooneyRivlinModel, "equations/E12")
CellML.VariableSetAsKnown(MooneyRivlinModel, "equations/E13")
CellML.VariableSetAsKnown(MooneyRivlinModel, "equations/E22")
CellML.VariableSetAsKnown(MooneyRivlinModel, "equations/E23")
CellML.VariableSetAsKnown(MooneyRivlinModel, "equations/E33")
# and variables to get from the CellML 
CellML.VariableSetAsWanted(MooneyRivlinModel, "equations/Tdev11")
CellML.VariableSetAsWanted(MooneyRivlinModel, "equations/Tdev12")
CellML.VariableSetAsWanted(MooneyRivlinModel, "equations/Tdev13")
CellML.VariableSetAsWanted(MooneyRivlinModel, "equations/Tdev22")
CellML.VariableSetAsWanted(MooneyRivlinModel, "equations/Tdev23")
CellML.VariableSetAsWanted(MooneyRivlinModel, "equations/Tdev33")

Having flagged the variables we require from the CellML model, we can map them to fields, more specifically, the components of field variables. First we map the strain tensor from the finite elasticity model to variables in the CellML model:

1
2
3
4
5
6
7
#Map the strain components
CellML.CreateFieldToCellMLMap(dependentField,CMISS.FieldVariableTypes.U1,1, CMISS.FieldParameterSetTypes.VALUES,MooneyRivlinModel,"equations/E11", CMISS.FieldParameterSetTypes.VALUES)
CellML.CreateFieldToCellMLMap(dependentField,CMISS.FieldVariableTypes.U1,2, CMISS.FieldParameterSetTypes.VALUES,MooneyRivlinModel,"equations/E12", CMISS.FieldParameterSetTypes.VALUES)
CellML.CreateFieldToCellMLMap(dependentField,CMISS.FieldVariableTypes.U1,3, CMISS.FieldParameterSetTypes.VALUES,MooneyRivlinModel,"equations/E13", CMISS.FieldParameterSetTypes.VALUES)
CellML.CreateFieldToCellMLMap(dependentField,CMISS.FieldVariableTypes.U1,4, CMISS.FieldParameterSetTypes.VALUES,MooneyRivlinModel,"equations/E22", CMISS.FieldParameterSetTypes.VALUES)
CellML.CreateFieldToCellMLMap(dependentField,CMISS.FieldVariableTypes.U1,5, CMISS.FieldParameterSetTypes.VALUES,MooneyRivlinModel,"equations/E23", CMISS.FieldParameterSetTypes.VALUES)
CellML.CreateFieldToCellMLMap(dependentField,CMISS.FieldVariableTypes.U1,6, CMISS.FieldParameterSetTypes.VALUES,MooneyRivlinModel,"equations/E33", CMISS.FieldParameterSetTypes.VALUES)

In the CellML model, each component of the strain tensor has a variable: equations/E11, equations/E12, etc.; following the CellML variable name convention. In the finite elasticity model, the strain tensor is found in the dependent field as the U1 field variable. These mappings are defined using the CellML.CreateFieldToCellMLMap() method as the value of these variables in the CellML model is known in the finite elasticity model and will be set by OpenCMISS when evaluating the CellML model.

In a similar manner we define equivalent mappings for the stress tensor (the U2 field variable in the dependent field):

1
2
3
4
5
6
7
#Map the stress components
CellML.CreateCellMLToFieldMap(MooneyRivlinModel,"equations/Tdev11", CMISS.FieldParameterSetTypes.VALUES,dependentField,CMISS.FieldVariableTypes.U2,1,CMISS.FieldParameterSetTypes.VALUES)
CellML.CreateCellMLToFieldMap(MooneyRivlinModel,"equations/Tdev12", CMISS.FieldParameterSetTypes.VALUES,dependentField,CMISS.FieldVariableTypes.U2,2,CMISS.FieldParameterSetTypes.VALUES)
CellML.CreateCellMLToFieldMap(MooneyRivlinModel,"equations/Tdev13", CMISS.FieldParameterSetTypes.VALUES,dependentField,CMISS.FieldVariableTypes.U2,3,CMISS.FieldParameterSetTypes.VALUES)
CellML.CreateCellMLToFieldMap(MooneyRivlinModel,"equations/Tdev22", CMISS.FieldParameterSetTypes.VALUES,dependentField,CMISS.FieldVariableTypes.U2,4,CMISS.FieldParameterSetTypes.VALUES)
CellML.CreateCellMLToFieldMap(MooneyRivlinModel,"equations/Tdev23", CMISS.FieldParameterSetTypes.VALUES,dependentField,CMISS.FieldVariableTypes.U2,5,CMISS.FieldParameterSetTypes.VALUES)
CellML.CreateCellMLToFieldMap(MooneyRivlinModel,"equations/Tdev33", CMISS.FieldParameterSetTypes.VALUES,dependentField,CMISS.FieldVariableTypes.U2,6,CMISS.FieldParameterSetTypes.VALUES)

In this case we use the CellML.CreateCellMLToFieldMap() method as it is the CellML model which defines the calculation of stress for a given strain state (i.e., the stress tensor components are wanted variables from the CellML model). Defining these field mappings informs OpenCMISS that when the CellML model is evaluated known variables should have their values updated to the current state of the corresponding field variable components and following an evaluation the field variable components mapped to wanted variables should be updated to reflect the newly computed value of the variables in the CellML model.

Now that we have imported the CellML model that we wish to use in our simulation and appropriately flagged the relevant variables, we can finish the creation of our CellML environment:

1
CellML.CreateFinish()

Finishing the CellML environment creation will now trigger the OpenCMISS to instantiate the CellML model(s) in that environment as executable code. No further changes to that code are possible.

We now need to define the CellML models field for our finite elasticity model. First (lines 2-4 below) we create a default OpenCMISS field and set it as the CellML environment’s models field. We then iterate over all elements in our finite element model and set each Gauss point in all elements to be associated with the Mooney-Rivlin model we imported into our CellML environment above.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
#Create the CellML models field
CellMLModelsField = CMISS.Field()
CellML.ModelsFieldCreateStart(CellMLModelsFieldUserNumber,CellMLModelsField)
CellML.ModelsFieldCreateFinish()

# Stress, strain fields are evaluated at Gauss Points and any fields used in the computation should also have values at Gauss points
# Iterate through each element and set it for the elements gauss points

xidiv = 1.0/(NumberOfGaussXi+1)
for elem in exregion.elements:
#Gauss point number counter    
    ctr = 0
#Assign model for each quadraturePoint:
    for xi in range(0,NumberOfGaussXi):
        xi1 = (1.0+xi)*xidiv
        for xj in range(0,NumberOfGaussXi):
            xi2 = (1.0+xj)*xidiv
            for xk in range(0,NumberOfGaussXi):
                xi3 = (1.0+xk)*xidiv
                ctr = ctr + 1
                CellMLModelsField.ParameterSetUpdateGaussPoint(CMISS.FieldVariableTypes.U,
                                                               CMISS.FieldParameterSetTypes.VALUES,
                                                               elem.number,
                                                               ctr,
                                                               1,
                                                               MooneyRivlinModel)

The CellML parameters field and CellML intermediate field are simply created with standard default OpenCMISS fields, as shown below.

1
2
3
4
5
6
7
8
9
#Create the CellML parameters field --- the strain field
CellMLParametersField = CMISS.Field()
CellML.ParametersFieldCreateStart(CellMLParametersFieldUserNumber,CellMLParametersField)
CellML.ParametersFieldCreateFinish()

#  Create the CellML intermediate field --- the stress field
CellMLIntermediateField = CMISS.Field()
CellML.IntermediateFieldCreateStart(CellMLIntermediateFieldUserNumber,CellMLIntermediateField)
CellML.IntermediateFieldCreateFinish()

The CellML environment is now set-up and ready to use with the single Mooney-Rivlin model that we have imported. When defining the finite elasticity problem in OpenCMISS it is important that the CellML subtype is specified to ensure that the CellML constitutive law is used:

1
2
3
4
5
6
7
#Define the problem
problem = CMISS.Problem()
problem.CreateStart(problemUserNumber)
problem.SpecificationSet(CMISS.ProblemClasses.ELASTICITY,
    CMISS.ProblemTypes.FINITE_ELASTICITY,
    CMISS.ProblemSubTypes.FINITE_ELASTICITY_CELLML)
problem.CreateFinish()

And then we need to define the numerical solver to use for the CellML constitutive law that we are using:

1
2
3
4
5
6
7
8
#Create the problem solver CellML equations
CellMLSolver = CMISS.Solver()
problem.CellMLEquationsCreateStart()
nonLinearSolver.NewtonCellMLSolverGet(CellMLSolver)
CellMLEquations = CMISS.CellMLEquations()
CellMLSolver.CellMLEquationsGet(CellMLEquations)
CellMLEquations.CellMLAdd(CellML)
problem.CellMLEquationsCreateFinish()

In line 2 we create a new default solver and then initialise it as the CellML solver from the general non-linear solver that we are using for the finite elasticity model (line 4). Lines 5 and 6 show the creation of the CellML equations and initialising them from the CellML solver. Finally, in line 7 we add the CellML environment containing our constitutive model to the CellML equations.

The CellML constitutive law is now defined as part of the general finite elasticity model and the simulation can be performed following assignment of the required boundary conditions, etc.

Results

_images/results.png

Figure: Results from running this pipe extension simulation. The gold lines show the original, undeformed, cylinder geometry. The coloured lines show the deformed geometry, with the colour varying to show the difference in strain through the wall of the cylinder. The cones represent the three normalised principal strains at material points throughout the tissue volume (red for compression and blue for extension).

Fluid Mechanics: Navier-Stokes: 1D-0D Visible Human Example

Section author: David Ladd

Introduction

Computational models can be used to gain mechanical insight into blood flowing in vessels under healthy and pathophysiological states. However, constructing full subject-specific CFD models of the entire arterial and/or venous vasculature is currently considered impractical, owing to: (1) the time and resources required to identify, segment, and constrain a model of the billions of vessels in a human body and (2) the computational cost such a model would incur.

As a result, the systemic context of a modeled system is often incorporated into definable models through boundary conditions and/or coupled models. These multiscale models restrict higher dimensional (3D/1D) models to a region of interest and rely on simpler, lower dimensional (0D), models to approximate physical behaviour outside the higher-dimensional domain. This involves coupling together of dependent fields (ie pressure and velocity/flow), material fields (eg fluid viscosity and wall compliance), and geometric fields (eg vessel diameter) at the interfaces between 3D, 1D, and/or 0D model domains.

In the following example, a 1D network of 24 major arteries has been constructed from the male Visible Human dataset, as shown below. The resulting mesh contains 135 nodes on 67 quadratic line elements, as shown below.

_images/1DVessels.png _images/1DVHMesh.png

Figure: 1D arteries and mesh.

Simple 0D/lumped parameter RCR Windkessel models are coupled to the 1D model at the outlet boundaries. These models use the fluid-electric analog to provide an basic approximation of resistance to flow due to perfusing vascular beds (the R1 and R2 components) and downstream vessel compliance (the C component).

_images/RCR.png

Figure: The 3-element RCR Windkessel model

Methods

A flow waveform from published data is applied at the aortic root of the model and interpolated in time from tabulated data using cubic splines. Outlet boundary conditions are provided by the 0D solution.

OpenCMISS/CellML field mapping capabilities allow for the coupling of the 1D (OpenCMISS) and 0D (CellML) solvers:

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
#================================================================================================================================
#  CellML Model Maps
#================================================================================================================================

if (RCRBoundaries):

    #----------------------------------------------------------------------------------------------------------------------------
    # Description
    #----------------------------------------------------------------------------------------------------------------------------
    # A CellML OD model is used to provide the impedance from the downstream vascular bed beyond the termination
    # point of the 1D model. This is iteratively coupled with the the 1D solver. In the case of a simple resistance
    # model, P=RQ, which is analogous to Ohm's law: V=IR. A variable map copies the guess for the FlowRate, Q at 
    # the boundary from the OpenCMISS Dependent Field to the CellML equation, which then returns presssure, P.
    # The initial guess value for Q is taken from the previous time step or is 0 for t=0. In OpenCMISS this P value is 
    # then used to compute a new Area value based on the P-A relationship and the Riemann variable W_2, which gives a
    # new value for Q until the values for Q and P converge within tolerance of the previous value.
    #----------------------------------------------------------------------------------------------------------------------------

    qCellMLComponent = 1
    pCellMLComponent = 2

    # Create the CellML environment
    CellML = CMISS.CellML()
    CellML.CreateStart(CellMLUserNumber,Region)
    # Number of CellML models
    CellMLModelIndex = [0]*(numberOfTerminalNodes+1)

    # Windkessel Model
    for terminalIdx in range (1,numberOfTerminalNodes+1):
        nodeIdx = coupledNodeNumber[terminalIdx]
        # Veins output
        if (nodeIdx%2 == 0):
            nodeDomain = Decomposition.NodeDomainGet(nodeIdx,meshComponentNumberSpace)
            if nodeDomain == computationalNodeNumber:
                CellMLModelIndex[terminalIdx] = CellML.ModelImport("./Input/CellMLModels/"+str(terminalIdx)+"/ModelHeart.cellml")
                # known (to OpenCMISS) variables
                CellML.VariableSetAsKnown(CellMLModelIndex[terminalIdx],"Heart/Qi")
                CellML.VariableSetAsKnown(CellMLModelIndex[terminalIdx],"Heart/Po")
                # to get from the CellML side 
                CellML.VariableSetAsWanted(CellMLModelIndex[terminalIdx],"Heart/Qo")
                CellML.VariableSetAsWanted(CellMLModelIndex[terminalIdx],"Heart/Pi")
        # Arteries outlet
        else:
            print('reading model: ' + "./Input/CellMLModels/"+str(terminalIdx)+"/ModelRCR.cellml")
            nodeDomain = Decomposition.NodeDomainGet(nodeIdx,meshComponentNumberSpace)
            if nodeDomain == computationalNodeNumber:
                CellMLModelIndex[terminalIdx] = CellML.ModelImport("./Input/CellMLModels/"+str(terminalIdx)+"/ModelRCR.cellml")
                # known (to OpenCMISS) variables
                CellML.VariableSetAsKnown(CellMLModelIndex[terminalIdx],"Circuit/Qin")
                # to get from the CellML side 
                CellML.VariableSetAsWanted(CellMLModelIndex[terminalIdx],"Circuit/Pout")
    CellML.CreateFinish()

    # Start the creation of CellML <--> OpenCMISS field maps
    CellML.FieldMapsCreateStart()
    
    # ModelIndex
    for terminalIdx in range (1,numberOfTerminalNodes+1):
        nodeIdx = coupledNodeNumber[terminalIdx]
        # Veins output
        if (nodeIdx%2 == 0):
            nodeDomain = Decomposition.NodeDomainGet(nodeIdx,meshComponentNumberSpace)
            if nodeDomain == computationalNodeNumber:
                # Now we can set up the field variable component <--> CellML model variable mappings.
                # Map the OpenCMISS boundary flow rate values --> CellML
                # Q is component 1 of the DependentField
                CellML.CreateFieldToCellMLMap(DependentFieldNavierStokes,CMISS.FieldVariableTypes.U,1,
                 CMISS.FieldParameterSetTypes.VALUES,CellMLModelIndex[terminalIdx],"Heart/Qi",CMISS.FieldParameterSetTypes.VALUES)
                CellML.CreateFieldToCellMLMap(DependentFieldNavierStokes,CMISS.FieldVariableTypes.U2,1,
                 CMISS.FieldParameterSetTypes.VALUES,CellMLModelIndex[terminalIdx],"Heart/Po",CMISS.FieldParameterSetTypes.VALUES)
                # Map the returned pressure values from CellML --> CMISS
                # pCellML is component 2 of the Dependent field U1 variable
                CellML.CreateCellMLToFieldMap(CellMLModelIndex[terminalIdx],"Heart/Qo",CMISS.FieldParameterSetTypes.VALUES,
                 DependentFieldNavierStokes,CMISS.FieldVariableTypes.U1,qCellMLComponent,CMISS.FieldParameterSetTypes.VALUES)
                CellML.CreateCellMLToFieldMap(CellMLModelIndex[terminalIdx],"Heart/Pi",CMISS.FieldParameterSetTypes.VALUES,
                 DependentFieldNavierStokes,CMISS.FieldVariableTypes.U1,pCellMLComponent,CMISS.FieldParameterSetTypes.VALUES)
        # Arteries output
        else:
            nodeDomain = Decomposition.NodeDomainGet(nodeIdx,meshComponentNumberSpace)
            if nodeDomain == computationalNodeNumber:
                # Now we can set up the field variable component <--> CellML model variable mappings.
                # Map the OpenCMISS boundary flow rate values --> CellML
                # Q is component 1 of the DependentField
                CellML.CreateFieldToCellMLMap(DependentFieldNavierStokes,CMISS.FieldVariableTypes.U,1,
                                              CMISS.FieldParameterSetTypes.VALUES,CellMLModelIndex[terminalIdx],"Circuit/Qin",CMISS.FieldParameterSetTypes.VALUES)
                # Map the returned pressure values from CellML --> CMISS
                # pCellML is component 1 of the Dependent field U1 variable
                CellML.CreateCellMLToFieldMap(CellMLModelIndex[terminalIdx],"Circuit/Pout",CMISS.FieldParameterSetTypes.VALUES,
                                              DependentFieldNavierStokes,CMISS.FieldVariableTypes.U1,pCellMLComponent,CMISS.FieldParameterSetTypes.VALUES)

    # Finish the creation of CellML <--> OpenCMISS field maps
    CellML.FieldMapsCreateFinish()

    CellMLModelsField = CMISS.Field()
    CellML.ModelsFieldCreateStart(CellMLModelsFieldUserNumber,CellMLModelsField)
    CellML.ModelsFieldCreateFinish()
    
    # Set the models field at boundary nodes
    for terminalIdx in range (1,numberOfTerminalNodes+1):
        nodeIdx = coupledNodeNumber[terminalIdx]
        nodeDomain = Decomposition.NodeDomainGet(nodeIdx,meshComponentNumberSpace)
        if nodeDomain == computationalNodeNumber:
            CellMLModelsField.ParameterSetUpdateNode(CMISS.FieldVariableTypes.U,CMISS.FieldParameterSetTypes.VALUES,
             versionIdx,derivIdx,nodeIdx,1,CellMLModelIndex[terminalIdx])

    CellMLStateField = CMISS.Field()
    CellML.StateFieldCreateStart(CellMLStateFieldUserNumber,CellMLStateField)
    CellML.StateFieldCreateFinish()

    CellMLParametersField = CMISS.Field()
    CellML.ParametersFieldCreateStart(CellMLParametersFieldUserNumber,CellMLParametersField)
    CellML.ParametersFieldCreateFinish()

    CellMLIntermediateField = CMISS.Field()
    CellML.IntermediateFieldCreateStart(CellMLIntermediateFieldUserNumber,CellMLIntermediateField)
    CellML.IntermediateFieldCreateFinish()

    # Finish the parameter update
    DependentFieldNavierStokes.ParameterSetUpdateStart(CMISS.FieldVariableTypes.U1,CMISS.FieldParameterSetTypes.VALUES)
    DependentFieldNavierStokes.ParameterSetUpdateFinish(CMISS.FieldVariableTypes.U1,CMISS.FieldParameterSetTypes.VALUES)

Snippet: OpenCMISS/CellML field mappings

Flow rate (Q) from the 1D Navier-Stokes/Characteristic solver provides the forcing term for the 0D ODE circuit solver. Pressure (P) is returned from the CellML solver to provide constraints on the Riemann invariants of the 1D system, which translate to area boundary conditions for the 1D FEM solver. At each timestep, the 1D and 0D systems are iteratively coupled until the boundary values converge within a user-specified tolerance at the 1D-0D interfaces. This procedure is outlined in the figure below.

_images/1D0DSolver.png

Figure: Overview of the coupled 1D-0D solution process.

Results

_images/1D0DVHFlowrates.png

Figure: Flow rates from the 1D-0D solution. Vessels shown at peak systole

_images/1D0DVHPressures.png

Figure: Pressure from the 1D-0D solution. Vessels shown at peak systole

Monodomain on a 2D domain with Noble 98 cellular model

Section author: Chris Bradley (c.bradley@auckland.ac.nz)

This example demonstrates a monodomain solution in a 2D domain using a Noble 98 cellular model. The Noble 98 model is used via CellML.

This is a Python example. In order to exectute this example, you will need to have the OpenCMISS Python bindings available. This Cmgui command file can be used to visualise the simulation output.

The example is encoded in Monodomain2DSquare.py. Below we describe some of the key aspects of this example.

Example description

This example starts, like all python examples with importing modules. In addition to sys, os and the CMISS module, math is imported as we will be using the math.sqrt function later on.

1
2
3
4
5
6
import sys, os, math
# Make sure $OPENCMISS_ROOT/cm/bindings/python is first in our PYTHONPATH.
sys.path.insert(1, os.path.join((os.environ['OPENCMISS_ROOT'],'cm','bindings','python')))

# Intialise OpenCMISS
from opencmiss import CMISS

The exact setup of the problem is controlled by a number of parameters in the python script, as shown below. The example solves the monodomain equations on a 2D domain of dimensions width by height. It is discretised with bilinear Lagrange finite elements with numberOfXElements in the X direction and numberOfYElements in the Y direction. The material parameters for the domain are controlled with the Am parameter for the cell membrane area, Cm for the membrane capacitance and conductivity for continuum conducticity. To start the simulation a stimulation is applied to the left half of the bottom row of nodes. This stimulation lasts from time zero until time stimStop. The magnitude of the stimulus is given by stimValue parameter. After time stimStop the stimulation is turned of and the simulation is continued until time timeStop. The time step for the spatial PDE problem is given by the pdeTimeStep parameter and the time step for the ODE integration is given by the odeTimeStep parameter. The final parameter, outputFrequency, controls how many time steps pass before the solution is output to file.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
# 2D domain size
height = 1.0
width = 1.0
numberOfXElements = 25
numberOfYElements = 25

# Materials parameters
Am = 193.6
Cm = 0.014651
conductivity = 0.1

# Simulation parameters
stimValue = 100.0
stimStop = 0.1
timeStop = 1.5
odeTimeStep = 0.00001
pdeTimeStep = 0.001
outputFrequency = 1

This example is designed to by run in parallel using MPI. OpenCMISS provides access to useful information about the parallel runtime environment via the calls below. numberOfComputationalNodes gives the total number of computational nodes the example is being run on. These nodes are numbered from 0 to (numberOfComputationalNodes - 1). computationalNodeNumber gives the node number that this particular process is running on.

1
2
3
# Get the number of computational nodes and this computational node number
numberOfComputationalNodes = CMISS.ComputationalNumberOfNodesGet()
computationalNodeNumber = CMISS.ComputationalNodeNumberGet()

The first step in our example will be to initialise OpenCMISS an to set up a region and coordinate system for our 2D mesh. This is achieved with the following calls

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
# Create a 2D rectangular cartesian coordinate system
coordinateSystem = CMISS.CoordinateSystem()
coordinateSystem.CreateStart(coordinateSystemUserNumber)
coordinateSystem.DimensionSet(2)
coordinateSystem.CreateFinish()

# Create a region and assign the coordinate system to the region
region = CMISS.Region()
region.CreateStart(regionUserNumber,CMISS.WorldRegion)
region.LabelSet("Region")
region.coordinateSystem = coordinateSystem
region.CreateFinish()

In this example we will use the generated mesh capabilities of OpenCMISS to generate our 2D mesh. The first thing we need to do is create a basis function for the mesh. We will define a bilinear Lagrange basis.

1
2
3
4
5
6
7
8
# Define a bilinear Lagrange basis
basis = CMISS.Basis()
basis.CreateStart(basisUserNumber)
basis.type = CMISS.BasisTypes.LAGRANGE_HERMITE_TP
basis.numberOfXi = 2
basis.interpolationXi = [CMISS.BasisInterpolationSpecifications.LINEAR_LAGRANGE]*2
basis.quadratureNumberOfGaussXi = [3]*2
basis.CreateFinish()

We can now create a generated mesh. We will create a regular mesh of size width x height and divide the mesh into numberOfXElements in the X direction and numberOfYElements in the Y direction. Note that when we finish generating the mesh we have a mesh object returned to us. This mesh object is just the same as if we had manually created the regular mesh.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
# Create a generated mesh
generatedMesh = CMISS.GeneratedMesh()
generatedMesh.CreateStart(generatedMeshUserNumber,region)
generatedMesh.type = CMISS.GeneratedMeshTypes.REGULAR
generatedMesh.basis = [basis]
generatedMesh.extent = [width,height]
generatedMesh.numberOfElements = [numberOfXElements,numberOfYElements]

mesh = CMISS.Mesh()
generatedMesh.CreateFinish(meshUserNumber,mesh)

Once the mesh has been created we can decompose it into a number of domains in order to allow for parallelism. We choose the options to let OpenCMISS calculate the best way to break up the mesh. We also set the number of domains to be equal to the number of computational nodes this example is running on.

1
2
3
4
5
6
# Create a decomposition for the mesh
decomposition = CMISS.Decomposition()
decomposition.CreateStart(decompositionUserNumber,mesh)
decomposition.type = CMISS.DecompositionTypes.CALCULATED
decomposition.numberOfDomains = numberOfComputationalNodes
decomposition.CreateFinish()

Now that the mesh has been decomposed we are in a position to create fields. The first field we need to create is the geometry field. Here we create a field and set the field’s mesh decomposition to the decomposed mesh that we have just created. We can choose exact how each component of the field is interpolated by setting component mesh component to be the mesh components that we created for the mesh. For this example we only have one mesh component. Once we have finished creating the field we can change the field DOFs to give us our geometry. Since this mesh has been generated we can use the generated mesh object to calculate the geometric parameters of the regular mesh.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
# Create a field for the geometry
geometricField = CMISS.Field()
geometricField.CreateStart(geometricFieldUserNumber, region)
geometricField.meshDecomposition = decomposition
geometricField.TypeSet(CMISS.FieldTypes.GEOMETRIC)
geometricField.VariableLabelSet(CMISS.FieldVariableTypes.U, "coordinates")
geometricField.ComponentMeshComponentSet(CMISS.FieldVariableTypes.U, 1, linearMeshComponentNumber)
geometricField.ComponentMeshComponentSet(CMISS.FieldVariableTypes.U, 2, linearMeshComponentNumber)
geometricField.CreateFinish()

# Set geometry from the generated mesh
generatedMesh.GeometricParametersCalculate(geometricField)

We are now in a position to define the type of physics that we wish to solve. This is done by creating an equations set which is a contianer object for all the parameters we need to describe the physics. Here we create an bioelectrics class equation of type monodomain.

1
2
3
4
5
6
7
8
9
# Create the equations_set
equationsSetField = CMISS.Field()
equationsSet = CMISS.EquationsSet()
equationsSet.CreateStart(equationsSetUserNumber, region, geometricField,
        CMISS.EquationsSetClasses.BIOELECTRICS,
        CMISS.EquationsSetTypes.MONODOMAIN_EQUATION,
        CMISS.EquationsSetSubtypes.NONE,
        equationsSetFieldUserNumber, equationsSetField)
equationsSet.CreateFinish()

Next we create the fields for the equations set. For the monodomain equation we need a dependent field (our solution) and a materials field (to descibe the physical material constant. Here we do not define a field before the create starts and so we let OpenCMISS create an appropriate dependent and materials field for the monodomain equations being described. Once the fields have been created we can set the field DOF values. In this example the domain is isotropic and so we set the four components of the materials field to be the Am, Cm and conductivity that we defined in the parameters section.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
# Create the dependent Field
dependentField = CMISS.Field()
equationsSet.DependentCreateStart(dependentFieldUserNumber, dependentField)
equationsSet.DependentCreateFinish()

# Create the materials Field
materialsField = CMISS.Field()
equationsSet.MaterialsCreateStart(materialsFieldUserNumber, materialsField)
equationsSet.MaterialsCreateFinish()

# Set the materials values
# Set Am
materialsField.ComponentValuesInitialise(CMISS.FieldVariableTypes.U,CMISS.FieldParameterSetTypes.VALUES,1,Am)
# Set Cm
materialsField.ComponentValuesInitialise(CMISS.FieldVariableTypes.U,CMISS.FieldParameterSetTypes.VALUES,2,Cm)
# Set conductivity
materialsField.ComponentValuesInitialise(CMISS.FieldVariableTypes.U,CMISS.FieldParameterSetTypes.VALUES,3,conductivity)
materialsField.ComponentValuesInitialise(CMISS.FieldVariableTypes.U,CMISS.FieldParameterSetTypes.VALUES,4,conductivity)
CellML

As mentioned above, this example uses the Noble 98 guinea-pig electrophysiology model defined in CellML. The guinea-pig ventricular cell model, originally developed by Noble et al in 1991, has been greatly extended to include accumulation and depletion of calcium in a diadic space between the sarcolemma and the sarcoplasmic reticulum where, according to contempory understanding, the majority of calcium-induced calcium release is triggered. A schematic of the model is shown below

_images/noble_1998a.png

Figure: Schematic of the Noble 98 guinea-pig ventricular cell model.

In this section we briefly describe how a CellML cellular model is used for a monodomain simulation. Further details and general description can be found in the OpenCMISS CellML Examples. First the CellML environment must be created and our Nobel 98 model imported into the environment:

1
2
3
4
5
# Create the CellML environment
cellML = CMISS.CellML()
cellML.CreateStart(cellMLUserNumber, region)
# Import a Nobel 98 cell model from a file
noble98Model = cellML.ModelImport("n98.xml")

With the CellML model imported, we are now able to flag the variables from the model that are known (by the monodomain pde model) and those variables that are wanted (by the monodomain pde model) from the CellML model. These are the variables from the CellML model that we will want to associate with fields in the OpenCMISS (pde continuum) model.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
# Now we have imported the model we are able to specify which variables from the model we want to set from openCMISS
cellML.VariableSetAsKnown(noble98Model, "fast_sodium_current/g_Na")
cellML.VariableSetAsKnown(noble98Model, "membrane/IStim")
# and variables to get from the CellML 
cellML.VariableSetAsWanted(noble98Model, "membrane/i_K1")
cellML.VariableSetAsWanted(noble98Model, "membrane/i_to")
cellML.VariableSetAsWanted(noble98Model, "membrane/i_K")
cellML.VariableSetAsWanted(noble98Model, "membrane/i_K_ATP")
cellML.VariableSetAsWanted(noble98Model, "membrane/i_Ca_L_K")
cellML.VariableSetAsWanted(noble98Model, "membrane/i_b_K")
cellML.VariableSetAsWanted(noble98Model, "membrane/i_NaK")
cellML.VariableSetAsWanted(noble98Model, "membrane/i_Na")
cellML.VariableSetAsWanted(noble98Model, "membrane/i_b_Na")
cellML.VariableSetAsWanted(noble98Model, "membrane/i_Ca_L_Na")
cellML.VariableSetAsWanted(noble98Model, "membrane/i_NaCa")

The transmembrane voltage is a variable that is used by both the continuum pde model and the CellML model. The variable is thus mapped as both known and wanted. This means that the continuum value of V is copied to the CellML model before evaluation and the CellML value of V is copied from the CellML model to the continuum pde model once the CellML model has finished being integrated. The other variables that we required to be flagged are the stimulus current and the transmembrane currents. The stimulus current is flagged as known as it will be controled by the pde model. The

Having flagged the variables we require from the CellML model, we can map them to fields, more specifically, the components of field variables. We map the transmembrane voltage to variables in the CellML model. Here Vm is needs to be passed from the OpenCMISS fields to the CellML fields before the CellML model is evaluated and then passed out from CellML to OpenCMISS after it is evaluated.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
# Start the creation of CellML <--> OpenCMISS field maps
cellML.FieldMapsCreateStart()
#Now we can set up the field variable component <--> CellML model variable mappings.
#Map Vm
cellML.CreateFieldToCellMLMap(dependentField,CMISS.FieldVariableTypes.U,1, CMISS.FieldParameterSetTypes.VALUES,noble98Model,"membrane/V", CMISS.FieldParameterSetTypes.VALUES)
cellML.CreateCellMLToFieldMap(noble98Model,"membrane/V", CMISS.FieldParameterSetTypes.VALUES,dependentField,CMISS.FieldVariableTypes.U,1,CMISS.FieldParameterSetTypes.VALUES)

#Finish the creation of CellML <--> OpenCMISS field maps
cellML.FieldMapsCreateFinish()

# Set the initial Vm values
dependentField.ComponentValuesInitialise(CMISS.FieldVariableTypes.U, CMISS.FieldParameterSetTypes.VALUES, 1,-92.5)

Note that we have initialised the OpenCMISS field for the transmembrane voltage to be the resting transmembrane voltage. We have set the wanted flag on the transmembrane current variables in the CellML model. These variables will be stored in the intermediates field that we create below. As these variables do not affect any other OpenCMISS field we do not need to map them.

We now need to create the CellML fields. We do this by creating the models field

1
2
3
4
#Create the CellML models field
cellMLModelsField = CMISS.Field()
cellML.ModelsFieldCreateStart(cellMLModelsFieldUserNumber, cellMLModelsField)
cellML.ModelsFieldCreateFinish()

the state field

1
2
3
4
#Create the CellML state field 
cellMLStateField = CMISS.Field()
cellML.StateFieldCreateStart(cellMLStateFieldUserNumber, cellMLStateField)
cellML.StateFieldCreateFinish()

and the parameters and intermediates field

1
2
3
4
5
6
7
8
9
#Create the CellML parameters field 
cellMLParametersField = CMISS.Field()
cellML.ParametersFieldCreateStart(cellMLParametersFieldUserNumber, cellMLParametersField)
cellML.ParametersFieldCreateFinish()

#  Create the CellML intermediate field 
cellMLIntermediateField = CMISS.Field()
cellML.IntermediateFieldCreateStart(cellMLIntermediateFieldUserNumber, cellMLIntermediateField)
cellML.IntermediateFieldCreateFinish()

Results

_images/start.png

Figure: Transmembrane voltage field immediately after the start of the simulation. A stimulus current has been applied to the first half of the bottom row of nodes.

_images/normal.png

Figure: Transmembrane voltage field after a fixed period of the simulation with a uniform gNa distribution at its normal value. The transmembrane voltage varies from -95 mV (blue) to +50 mV (red).

_images/gNaDistribution.png

Figure: Plot of gNa with a radial variation. gNa varies from its normal value of 3.855 x 10^-5 mS mm^-2 at the bottom left node (blue) to 300% of its normal value at the top right node (red).

_images/gNa.png

Figure: Transmembrane voltage field after a fixed period of the simulation with a radial gNa distribution from 100 to 300 % of its normal value (shown above). The transmembrane voltage varies from -95 mV (blue) to +50 mV (red).

Glossary

CellML environment
CellML Environment
The CellML Environment is the primary object in the OpenCMISS API which manages the use of CellML models in Iron. A single instance of a CellML environment is able to manage multiple CellML models (and multiple instances of the same model, if required).
CellML intermediate field
The OpenCMISS field which is used to contain the intermediate variables for the CellML models for a given CellML environment, these are wanted variables which are not state variables in the CellML model(s).
CellML models field
The field which associates CellML models within a given CellML environment to finite element models in OpenCMISS.
CellML parameters field
The OpenCMISS field which is used to define spatially varying known parameters in the CellML models for a given CellML environment.
CellML variable ID
CellML variable name
In order to uniquely identify variables from CellML models using a text string, OpenCMISS methods use identifiers of the form: component_name/variable_name. Where component_name is the value of the name attribute of the component in which the desired variable is located and the variable_name is the value of that variable’s name attribute. Variables in the CellML model which are connected can be addressed by any of the relevant CellML variable ID’s.
known
Known variables
CellML models used in OpenCMISS are required to be mathematically complete and suitably constrained. Often for such models to be useful in simulations, the definition of specific variables from the model(s) must be overridden in order to allow their value to be controlled during the simulation by OpenCMISS fields. Such variables are said to be known. When addressing variables from CellML models in the OpenCMISS API we use the standard CellML variable ID.
wanted
Wanted variables
In order for OpenCMISS to have access to retrive the value of variables from CellML models, those variables must be flagged to OpenCMISS. Such variables are said to be wanted. When addressing variables from CellML models in the OpenCMISS API we use the standard CellML variable ID.

Indices and tables