September 24, 2017, Sunday

# Generalizing to More Points

## Introduction

We are interested in examining how our model responds to an increase in the number of points of approximation. Logically, we should expect more "informative" differential equations, providing information as to the behavior of the system at an increased number of points for some spatial domain, but fundamentally, we should not be expecting any important new information through SVD analysis, provided we originally chose an appropriate number of points representative of the different regions within our spatial domain (for example, in our C5a production model, we chose three points: at the bacterial cell surface, in the bacterial capsule, and at the capsule-plasma interface).

We begin by examining a concrete case, in which only one additional approximation point is added within the spatial domain of our toy examples, to give a total of four discrete, evenly-spaced points in the spatial domain. We will then generalize how our model behaves for some $m$ points, including the effect as $m\rightarrow \infty$.

## Four Points

The addition of the fourth point does not at all change the means by which we generate our diffusion differential equations. We require six points, which we denote as 0, 1, 2, 3, 4, 5, of which points 1 - 4 reside in the spatial domain of interest and concentration information at points 0 and 5 are necessary for our approximation method. Again, we define h to be the distance between two adjacent points, so that if we denote H as the total distance in the spatial domain, $h={\frac {H}{3}}$. Thus, we have the following system:

${\frac {\partial \phi _{{1}}}{\partial t}}=D({\frac {(\phi _{{0}}-\phi _{{1}})+(\phi _{{2}}-\phi _{{1}})}{h^{{2}}}})$
${\frac {\partial \phi _{{2}}}{\partial t}}=D({\frac {(\phi _{{1}}-\phi _{{2}})+(\phi _{{3}}-\phi _{{2}})}{h^{{2}}}})$
${\frac {\partial \phi _{{3}}}{\partial t}}=D({\frac {(\phi _{{2}}-\phi _{{3}})+(\phi _{{4}}-\phi _{{3}})}{h^{{2}}}})$
${\frac {\partial \phi _{{4}}}{\partial t}}=D({\frac {(\phi _{{3}}-\phi _{{4}})+(\phi _{{5}}-\phi _{{4}})}{h^{{2}}}})$

Again, we examine Dirichlet-Dirchlet, Dirichlet-Neumann, and Neumann-Neumann boundary conditions at the edges of our spatial domain separately.

Dirichlet-Dirichlet Boundary Condition

If we define Dirchlet boundary conditions at both the left and right boundaries of our spatial domain, we have:

$\phi _{{0}}-\phi _{{1}}=0$ and
$\phi _{{5}}-\phi _{{4}}=0$,

eliminating these two fluxes from our diffusion matrix.

Additionally, we now have three equal and opposite fluxes of the form specified in our original three-point example, cancelling an additional three fluxes that need to be represented in our matrix. Thus, we have a total of four "reacting species" (four points) and three remaining diffusive fluxes. Therefore, we construct a $4\times 3$ diffusion matrix D:

${\textbf {D}}_{{Dirichlet-Dirchlet}}={\begin{pmatrix}1&0&0\\-1&1&0\\0&-1&1\\0&0&-1\end{pmatrix}}$
where
$\phi =(\phi _{{1}},\phi _{{2}},\phi _{{3}},\phi _{{4}})$
and
${\textbf {v}}=(D{\frac {\phi _{{2}}-\phi _{{1}}}{h^{{2}}}},D{\frac {\phi _{{3}}-\phi _{{2}}}{h^{{2}}}},D{\frac {\phi _{{4}}-\phi _{{3}}}{h^{{2}}}})^{{{\textbf {T}}}}$

Let $m=4$, corresponding to the number of rows of our matrix, and $n=3$, corresponding to the number of columns of our matrix. Using MATLAB's rank operation for $A$ gives a rank $r=3$.

The dimension of the left null space is $m-r=4-3=1$, which means that the left null space has one spanning vector. Through our SVD analysis, we see from the matrix $U$ that an orthonormal basis for the left null space is composed of the vector $U_{{4}}$ alone.

$U_{{4}}=({\frac {1}{2}},{\frac {1}{2}},{\frac {1}{2}},{\frac {1}{2}})$

This corresponds to a conserved quantity, where $\phi _{{1}}+\phi _{{2}}+\phi _{{3}}+\phi _{{4}}$ is time invariant. As expected, this result is analogous to our result for three approximation points, with the sum of the reacting species at the discrete points constant in time.

Meanwhile, the dimension of the right null space is $n-r=3-3=0$, meaning that the right null space has no spanning vectors. Biologically, there are no steady-state flux distributions, besides the trivial state in which the initial reactant concentrations are equal and no diffusive gradients form. Again, this result is functionally the same as that for three approximation points.

Dirichlet-Neumann Boundary Conditions

If we define a Dirchlet boundary condition at the left boundary of our spatial domain and a Neumann boundary condition at the right boundary, we have:

$\phi _{{0}}-\phi _{{1}}=0$ and
$D({\frac {\phi _{{5}}-\phi _{{4}}}{h^{{2}}}})=\beta$,

giving us an additional flux in our system compared to the Dirichlet-Dirichlet case.

Thus, we can construct a $4\times 4$ matrix:

${\textbf {D}}_{{Dirichlet-Neumann}}={\begin{pmatrix}1&0&0&0\\-1&1&0&0\\0&-1&1&0\\0&0&-1&1\end{pmatrix}}$
where
$\phi =(\phi _{{1}},\phi _{{2}},\phi _{{3}},\phi _{{4}})$
and
${\textbf {v}}=(D{\frac {\phi _{{2}}-\phi _{{1}}}{h^{{2}}}},D{\frac {\phi _{{3}}-\phi _{{2}}}{h^{{2}}}},D{\frac {\phi _{{4}}-\phi _{{3}}}{h^{{2}}}},\beta )^{{{\textbf {T}}}}$

Let $m=4$, corresponding to the number of rows of our matrix, and $n=4$, corresponding to the number of columns of our matrix. MATLAB's rank operation for $A$ gives a rank $r=4$.

The dimension of the left null space is $m-r=4-4=0$, which means that the left null space has no spanning vectors. Biologically, there are no conserved quantities, which is the same result that we obtained for the three-point situation.

Meanwhile, the dimension of the right null space is $n-r=4-4=0$, meaning that the right null space has no spanning vectors. Biologically, there are no steady-state flux distributions, besides the trivial state in which the initial reactant concentrations are equal and no diffusive gradients form. Again, this result is equivalent to that obtained for the three point situation.

Neumann-Neumann Boundary Conditions

If we define Neumann boundary conditions at both the left and right boundaries of our spatial domain, we have::

$D({\frac {\phi _{{0}}-\phi _{{1}}}{h^{{2}}}})=\beta _{{1}}$ and
$D({\frac {\phi _{{5}}-\phi _{{4}}}{h^{{2}}}})=\beta _{{2}}$,

with the subscripts 1 and 2 to distinguish the two potentially different fluxes at the left and right boundaries.

We have added yet another flux to the system relative to the Dirichlet-Neumann case. Thus, we can construct a $4\times 5$ matrix:

${\textbf {D}}_{{Neumann-Neumann}}={\begin{pmatrix}1&1&0&0&0\\0&-1&1&0&0\\0&0&-1&1&0\\0&0&0&-1&1\end{pmatrix}}$
where
$\phi =(\phi _{{1}},\phi _{{2}},\phi _{{3}},\phi _{{4}})$
and
${\textbf {v}}=(\beta _{{1}},D{\frac {\phi _{{2}}-\phi _{{1}}}{h^{{2}}}},D{\frac {\phi _{{3}}-\phi _{{2}}}{h^{{2}}}},D{\frac {\phi _{{4}}-\phi _{{3}}}{h^{{2}}}},\beta _{{2}})^{{{\textbf {T}}}}$

Let $m=4$, corresponding to the number of rows of our matrix, and $n=5$, corresponding to the number of columns of our matrix. MATLAB's rank operation for $A$ gives a rank $r=4$.

The dimension of the left null space is $m-r=4-4=0$, which means that the left null space has no spanning vectors. Biologically, there are no conserved quantities, which again matches our result for the approximation with three points.

Meanwhile, the dimension of the right null space is $n-r=5-4=1$, meaning that the right null space has one spanning vector. Through our SVD analysis, we see that an orthonormal basis is composed of the vector $V_{{5}}$ alone.

We reproduce this vector here using MATLAB's svd function:

$V_{{4}}=(-{\frac {1}{{\sqrt {5}}}},{\frac {1}{{\sqrt {5}}}},{\frac {1}{{\sqrt {5}}}},{\frac {1}{{\sqrt {5}}}},{\frac {1}{{\sqrt {5}}}})$

Remember that the Neumann boundary conditions were specified as +1 in the diffusive matrix, which corresponds to fluxes into the spatial domain. Biologically, our 'V' vector corresponds to the steady state distribution:

$-\beta _{{1}}=D{\frac {\phi _{{2}}-\phi _{{1}}}{h^{{2}}}}=D{\frac {\phi _{{3}}-\phi _{{2}}}{h^{{2}}}}=D{\frac {\phi _{{4}}-\phi _{{3}}}{h^{{2}}}}=\beta _{{2}}$

This makes biological sense, as it corresponds to the situation in which each point on the spatial domain is balanced by equal and opposite diffusion fluxes. Specifically, this steady-state case really corresponds to the reactant in question diffusing across the entirety of the spatial domain at a uniform rate. This result again corresponds to the one obtained for the three-point case.

## Infinite Points

As $m\rightarrow \infty$, we obtain an infinite matrix of both rows and columns, as for each approximating point added, our diffusive matrix must account a diffusive flux in to and out of one new point, generating an additional column, corresponding to an additional diffusive flux as well as an additional row, corresponding to an additional differential equation for an additional approximated point.

From our examples using $m=3$ and $m=4$, it is easy to see that our matrices follow a certain structure based on the number of points we choose to approximate at. As mentioned previously, we can consider the diffusion between points in our system as a network of reversible chemical reactions, where:

$\phi _{{i}}\leftrightarrow \phi _{{i+1}}$

for $i=0,2,...,m$ across our spatial domain and in agreement with Fick's First Law of diffusion. We know that stoichiometric matrices, read "down" the column, correspond to these chemical reactions, so that each column corresponding to a diffusive flux will have a '+1' followed by a '-1', corresponding to diffusion from point 'i' to point 'i + 1'. These entries will appear in the row to which the diffusion point corresponds to -- that is, as we have constructed our example matrices, the '+1' appears in row 'i', and the '-1' appears directly below it in row 'i + 1' for the diffusion flux from 'i' to 'i + 1'.

Of course, diffusion is by nature reversible, but as demonstrated when first constructing this matrix, the "flux" between two given points is equal and opposite when taken in reverse directions. Thus, adding an extra row to correspond to the reverse "reaction" of diffusion with a '-1' in row 'i' and a '+1' in row 'i + 1' is both biologically and mathematically redundant.

Dirichlet-Dirichlet Case

In the Dirchlet-Dirichlet boundary condition case, then, we obtain a matrix with the following general structure:

$D_{{Dirchlet-Dirichlet,n\rightarrow \infty }}={\begin{pmatrix}1&0&\cdots &0&0\\-1&1&\cdots &0&0\\0&-1&\cdots &0&0\\\vdots &\vdots &\ddots &\vdots &\vdots \\0&0&\cdots &1&0\\0&0&\cdots &-1&1\\0&0&\cdots &0&-1\\\end{pmatrix}}$

It is evident that this matrix is full rank, and thus of infinite rank as $m\rightarrow \infty$, as the columns corresponding to the diffusion reactions between separate points are linearly (and biologically) independent. Specifically, $r=n$, where n is the number of columns in the matrix, corresponding to the number of diffusion reactions (columns) described between points (reversible). In general, we know that $m=n+1$, as the diffusion reactions/columns represent the interactions between the points (which are now continuous). Thus, the total number of approximate points must be one greater than the total number of interactions between them.

For the infinite case, this rank implies that the dimension of the left null space is $m-n=1$: same as in our $m=3$ and $m=4$ cases. By inspection, it is evident that each entry in this left null space vector must have equal value. This is because each column contains equal but opposite entries corresponding to the diffusion from one point to another with "overlapping" columns, as the reactant diffuses from one point to another. Mathematically, this requires that the first and second entries be equal, the second and third entries be equal, and so on to produce a spanning vector in the left null space. In other words, because these spanning vectors can be regarded as a linear combination of nodes that sum to zero, and columns correspond to diffusive reactions between adjacent points, equal fluxes from point to point will give rise to a biological steady-state.

Meanwhile, there are no nontrivial spanning vectors in the null space, as we know that $n-r=0$ with a full rank matrix. Biologically, there are no steady-state flux distributions.

Dirichlet-Neumann Case

The Dirichlet-Neumann case differs from the Dirichlet-Dirichlet case only in that we introduce an additional "flux" between points two "outermost" points at one boundary of the spatial domain. As $m\rightarrow \infty$, this is effectively taken at a single point as traditionally understood. As seen in the examples for $m=3$ and $m=4$, this is mathematically represented through the addition of an extra column with some defined Neumann flux, say $\beta$, at the point for which the boundary condition is given. Since this column specifies some constant flux flowing in to (or out of) the spatial domain, it will only have a single entry and therefore be linearly independent of the other vectors found in the above section. Thus, this matrix is full rank -- one greater than the matrix described in the previous section. The number of rows $m$ remains the same, as we are taking the same "number" of infinite points in this continuous case, whereas the number of columns increases by one, as we are defining an additional flux at one boundary.

Therefore, there are no nontrivial spanning vectors in the left null space, as $m-r=0$, where $r$ has a value one greater here than in the Dirichlet-Dirichlet case. Biologically, there are no conserved quantities in this case.

Also, there are no nontrivial spanning vectors in the null space, as $n-r=0$ with a full rank matrix. Biologically, there are no steady-state flux distributions in this system.

Neumann-Neumann Case

The Neumann-Neumann case contains one more "flux" between the two "outermost" points at the second boundary of the spatial domain versus the Dirichlet-Neumann case, which with $m\rightarrow \infty$ is again essentially a boundary condition defined at a single point. Again, this is mathematically represented by adding an additional column with some defined Neumann flux at the point for which the other boundary condition is given. This column also has a single entry in the row that corresponds to the point at which the boundary condition is taken. This matrix is still full rank and equal to the rank of the Dirichlet-Neumann case, as no additional linearly independent rows were added (although a linearly independent column was) to a square matrix in the Dirichlet-Neumann case.

Again, there are no nontrivial spanning vectors in the left null space, as $m-r=0$ with no columns added from the previous case. Biologically, there are no conserved quantities in this case.

However, there is one spanning vector in the null space, since $n-r=1$ with an additional column corresponding with the Neumann flux at the other boundary of the spatial domain. The makeup of this particular vector depends on if we define flux in or out of the two boundaries. In general, though, since each row defines compounds such that each reactant point experiences two fluxes from opposite directions, the fluxes themselves must be equal but opposite to give rise to a steady-state flux distribution. Due to the overlapping nature of these fluxes, this implies that each flux in this spanning vector must have equal value. The signs on the interior of the spatial domain must be the same, as seen in the example matrices for $m=3$ and $m=4$, but at the boundaries, the sign of the fluxes depends on whether reactant is specified to flow in to or out of the boundary. Either way, the fluxes must balance such that the boundary points experience equal flux in to and out of the spatial domain.