Commit dadad563 authored by Félix Hartmann's avatar Félix Hartmann
Browse files

[bugfix] Fluxes at the boundary of the polar zone was erroneous.

parent 8b92cddb
......@@ -164,6 +164,9 @@ class CellFile(object):
(enl_id_conc >= self.enlargement_threshold)
& (div_id_conc < self.division_threshold))[0]
zone_3 = np.where(enl_id_conc < self.enlargement_threshold)[0]
zone_1.sort()
zone_2.sort()
zone_3.sort()
return zone_1, zone_2, zone_3
def get_cells_identities(self):
......@@ -240,7 +243,7 @@ class CellFile(object):
zone_2.append(cell.indices[i])
elif identity == "non growing":
zone_3.append(cell.indices[i])
return zone_1, zone_2, zone_3
return np.sort(zone_1), np.sort(zone_2), np.sort(zone_3)
def get_cell_concentrations(self, u, concentration_operator):
"""Return signal concentrations of the cells."""
......@@ -674,24 +677,30 @@ class CellFile(object):
polar_zone = np.arange(0, N)
polar_zone_lower = np.arange(0, N-1)
polar_zone_upper = np.arange(1, N)
# diagonal
D = 1 - 2*q * dt / self.D_real - d * dt - self.S * dt
D[0] += q * dt / self.D_real[0] # no diffusion at boundary
D[-1] += q * dt / self.D_real[-1] # no diffusion at boundary
# lower diagonal
Dlower = q * dt / self.D_real[1:]
if transport.carrier_orientation == "toward the xylem":
Dlower[polar_zone_lower] = \
(p + q) * dt / self.D_real[polar_zone_lower + 1]
np.fill_diagonal(M_i[1:, :-1], Dlower)
# upper diagonal
Dupper = q * dt / self.D_real[:-1]
# additional terms for polar transport
D[polar_zone] -= \
p * dt / self.D_real[polar_zone]
if transport.carrier_orientation == "toward the xylem":
Dlower[polar_zone_lower] += \
p * dt / self.D_real[polar_zone_lower + 1]
if N-1 in polar_zone: # no polar transport at boundary
D[-1] += p * dt / self.D_real[-1]
if transport.carrier_orientation == "toward the cambium":
Dupper[polar_zone_upper - 1] = \
(p + q) * dt / self.D_real[polar_zone_upper - 1]
np.fill_diagonal(M_i[:-1, 1:], Dupper)
# diagonal
D = 1 - 2*q * dt / self.D_real - d * dt - self.S * dt
D[polar_zone] = \
1 - (p + 2*q) * dt / self.D_real[polar_zone] \
- d * dt - self.S[polar_zone] * dt
Dupper[polar_zone_upper - 1] += \
p * dt / self.D_real[polar_zone_upper - 1]
if 0 in polar_zone: # no polar transport at boundary
D[0] += p * dt / self.D_real[0]
np.fill_diagonal(M_i, D)
np.fill_diagonal(M_i[1:, :-1], Dlower)
np.fill_diagonal(M_i[:-1, 1:], Dupper)
else:
RdX_flux_i = 1/self.D_real[self.flux_cells]
M_i = None
......
......@@ -301,7 +301,7 @@ def fill_cell_file_drawing(
# an arrow on one well-selected membrane
elif show_transport_polarity == "one":
zone_1, zone_2, zone_3 = cf.get_cell_zones_anatomy(time)
if zone_2:
if zone_2.size:
cell_index = zone_2[len(zone_2)//2]
else:
cell_index = nb_cells / 2
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment